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ABSTRACT 


Analyses  for  laminar  plates  are  reviewed  in  this  report. 
It  is  observed  that  while  approximate  global  stress  analyses 
have  been  performed  for  laminates,  little  has  been  accom¬ 
plished  in  the  area  of  strength  analysis.  The  work  presented 
here  is  concentrated  in  that  area,  i.e.,  mathematical  models 
to  consider  the  influence  of  flaws  on  failure  loads  asso¬ 
ciated  with  crack  propagation  and/or  delamination.  In  par¬ 
ticular,  an  approximate  three-dimensional  formulation  for 
laminates  based  on  the  variational  principle  of  minimum  com¬ 
plementary  potential  energy  is  presented  and  applied  to  the 
problem  of  a  three  layer  composite  plate  containing  a  through 
crack.  This  model  is  employed  to  study  the  influences  of 
the  geometric  and  material  characteristics  of  the  composite 
on  the  local  stress  field  in  the  vicinity  of  the  flaw. 

Results  of  an  experimental  program  to  determine  the  in¬ 
fluence  of  plate  thickness  on  critical  load  for  a  centrally 
cracked  tensile  specimen  are  reported.  Extension  of  the  ex¬ 
perimental  program  to  consideration  of  laminated  specimens  is 
currently  in  progress. 

A  stress  analysis  has  been  carried  out  for  another  prob¬ 
lem  of  particular  interest  in  the  study  of  failure  for  com¬ 
posite  systems.  The  effect  of  voids  in  adhesive  layers  on 
the  behavior  of  multi-layer  fibrous  composites  is  modeled  by 
the  problem  of  a  penny  shaped  flaw  contained  in  an  isotropic 


iii 


layer  joining  two  orthotropic  half  spaces.  Stresses  in  the 
vicinity  of  the  crack  edge  are  calculated  in  terms  of  the 
loading  and  the  properties  of  the  components  for  this  com¬ 
posite  system. 


iv 


Blank  Pages 


TABLE  OP  CONTENTS 


Section 

I.  INTRODUCTION 

II.  MODELING  OP  LAMINATED  PLATES 

2.1  Elastic  Properties  of  Composite  Materials 

2.2  Composite  Plates 

2.2.1  Sandwich  Plates 

2.2.2  Laminated  Composite  Plates 

2.2.3  N-Layer  Laminates 

2.2.4  Stress  Analysis  for  Laminates 
Containing  Flaws  -  A  Proposed 
Theory 

III.  AN  APPROXIMATE  THREE-DIMENSIONAL  LAMINAR  PLATE 
THEORY 

3*1  The  Variational  Principle  of  Minimum 
Complementary  Potential  Energy  for 
Composites 

3.2  An  Approximate  Three-Dimensional  Theory 
for  Laminates 

3.3  Application  of  the  Approximate  Theory 
to  a  Laminar  Plate  Containing  a  Through 
Crack 

3.3*1  Solution  Procedure  for  the  In- 
Plane  Variations  of  the  Stress 
Field 

3.3.2  Through-the-Thickness  Variation 
of  the  Near-Tip  Stress  Field 

3.3.3  Numerical  Results  and  Discussion 

3.3.4  Limitations  of  Present  Theory 
and  Future  Work 

IV.  EXPERIMENTAL  MEASUREMENT  OF  THE  EFFECT  OF  PLATE 
THICKNESS  ON  THE  STRAIN  ENERGY  RELEASE  RATE 

4.1  Calculation  of  Energy  Release  Rate  from 
Experimental  Results 

vi 


Page 

1 

5 

6 
8 
8 

10 

17 

31 


32 

34 


39 

52 


55 


68 

74 

79 

81 


82 


4.2  The  Experimental  Program  84 

4.3  Testing  Procedure  84 

4.4  Experimentally  Determined  Energy  Release  86 

Rates 

4.5  Comparison  of  Experimentally  Determined  86 

Energy  Release  Rates  with  the  Hartranft- 

Sih  Model 

V.  A  PENNY-SHAPED  CRACK  IN  AN  ISOTROPIC  LAYER  SAND-  94 
WICHED  BETWEEN  TWO  ORTHOTROPIC  HALF-SPACES 

5.1  Model  for  Fiber  Reinforced  Laminate  95 

5.2  Formulation  of  Governing  Equations  96 

5.3  Local  Stress  Field  106 

5.4  Numerical  Results  107 

5.5  Discussion  108 

APPENDIX  A  -  PARAMETERS  ARISING  IN  SOLUTION  FOR  IN-  110 

PLANE  STRESS  VARIATIONS 

APPENDIX  B  -  FUNCTION  f(£)  RELATED  TO  THE  SECOND  113 

EQUATION  (5.2.9) 

REFERENCES  118 

FIGURES  121 


vii 


LIST  OP  FIGURES 


Figure  3* 1.1  A  General  Two-Phase  Composite 

Figure  3.2.1  Laminated  Plate 

Figure  3.3*1  Three  Layer  Composite  Plate  Containing  a 

Through  Crack 

Figure  3 -  3 • 3 - 1  Dimensionless  Stress  Intensity  Factor  as  a 

Function  of  Normalized  Layer  Thickness  for 
Formulation  without  Boundary  Layer 

Figure  3 • 3 • 3 • 2  Dimensionless  Stress  Intensity  Factor  as  a 

Function  of  Normalized  Layer  Thickness  for 
Formulation  without  Boundary  Layer 

Figure  3. 3* 3* 3  Dimensionless  Stress  Intensity  Factor  as  a 

Function  of  Normalized  Layer  Thickness  for 
Formulation  without  Boundary  Layer 

Figure  3 -  3 • 3 • ^  Dimensionless  Stress  Intensity  Factor  as  a 

Function  of  Normalized  Layer  Thickness  for 
Formulation  without  Boundary  Layer 

Figure  3  -  3  -  3 • 5  Dimensionless  Stress  Intensity  Factor  at 

Various  Ratios  of  Material  Stiffness  for 
Formulation  including  Boundary  Layers 

Figure  3. 3. 3. 6  Dimensionless  Stress  Intensity  Factor  at 

Various  Ratios  of  Material  Stiffness  for 
Formulation  including  Boundary  Layers 

Figure  3- 3. 3* 7  Dimensionless  Stress  Intensity  Factor  at 

Various  Ratios  of  Material  Stiffness  for 
Formulation  including  Boundary  Layers 

Figure  3  -  3  -  3  -  8  Dimensionless  Stress  Intensity  Factor  at 

Various  Ratios  of  Material  Stiffness  for 
Formulation  including  Boundary  Layers 

Figure  3* 3* 3* 9  Dimensionless  Stress  Intensity  Factor  at 

Various  Ratios  of  Material  Stiffness  for 
Formulation  including  Boundary  Layers 

Figure  3. 3- 3. 10  Influence  of  Relative  Layer  Thickness  on 

Normalized  Stress  Intensity  Factor 

Figure  3. 3 -3. 11  Influence  of  Relative  Layer  Thickness  on 

Normalized  Stress  Intensity  Factor 

Figure  3.3.3.12  Influence  of  Relative  Layer  Thickness  on 

Normalized  Stress  Intensity  Factor 
viii 


Figure  3- 3 *3* 13  Influence  of  Relative  Layer  Thickness  on 

Normalized  Stress  Intensity  Factor 


Figure  3*3.3.14  Influence  of  Relative  Layer  Stiffness  and 


Thickness  on  Stress  Intensity  Factor  for 
Constant  Normalized  Plate  Thickness 

Figure 

4.1 

A  Central  Crack  in  an  Infinite  Plate 

Figure 

4.2a 

Plate  Containing  a  Central  Crack 

Figure 

4.2b 

Load-Displacement  Curve  for  Constant  Displace¬ 
ment  during  Crack  Extension 

Figure 

4.2c 

Load-Displacement  Curve  for  Constant  Load 
during  Crack  Extension 

Figure 

4.3 

Specimen  Geometry 

Figure 

4.4 

Testing  Details 

Figure 

4.5 

Testing  Details 

Figure 

4.6 

Load-Displacement  Curves 

Figure 

4.7 

Critical  Load  as  a  Function  of  Crack  Length 

Figure 

4.8 

Compliance  Curves  for  3/8"  Specimens 

Figure 

4.9 

Strain  Energy  Release  Rate  as  a  Function  of 
Normalized  Plate  Thickness 

Figure 

4.10 

Comparison  of  Predicted  and  Observed  Critical 
Loads  at  Various  Ratios  of  Plate  Thickness  to 
Crack  Length 

Figure 

5.1 

A  Penny-Shaped  Crack  in  an  Adhesive  Layer 
Between  Fiber  Reinforced  Layers 

Figure 

5.2 

Remote  Loading  for  the  Composite  Cylinder 
Containing  a  Penny-Shaped  Crack  in  the 

Adhesive  Layer 

Figure 

5-3 

Idealized  Model  for  the  Composite  Containing 
a  Penny-Shaped  Crack 

Figure 

5-4 

Variation  of  the  Normalized  Stress  Intensity 
Factor  with  Material  Properties  and  Relative 
Layer  Thickness 

Figure 

5-5 

Variation  of  the  Normalized  Stress  Intensity 
Factor  with  Material  Properties  and  Relative 
Layer  Thickness 

Figure  5*6 


Normalized  Stress  Intensity  Factor  as  a 
Function  of  Relative  Layer  Height 


SECTION  I 


INTRODUCTION 

A  research  effort  primarily  concerned  with  the  under¬ 
standing  and  prediction  of  failure  for  laminar  composites  is 
reported  here.  The  current  interest  in  fibrous  laminates 
for  structural  application  is  associated  with  the  high 
strength-to-weight  ratio  which  can  be  developed  in  laminates. 
These  laminates  are  generally  composed  of  layers  which  have 
been  reinforced  by  embedding  unidirectional  fibers.  The 
layers  are  adhered  to  each  other  such  that  the  fiber  direc¬ 
tion  varies  from  one  layer  to  the  next  in  a  previously  deter¬ 
mined  manner.  The  freedom  of  choice  for  fiber  orientation 
in  the  layers  of  the  composite  system  enables  the  development 
of  laminates  with  special  preferential  directional  properties 
for  particular  applications.  Because  of  this  characteristic 
of  fibrous  composites,  the  employment  of  these  systems  rather 
than  equivalent  homogeneous  members  will  be  clearly  advanta¬ 
geous  in  many  applications. 

The  drawback  of  using  laminates  is  that  they  are  sensitive 
to  a  number  of  manufacturing  imperfections  specific  to  com¬ 
posite  systems,  i.e.,  fiber  misalignment,  loading  imperfec¬ 
tions  between  the  fibers  and  matrix  material  within  an  indi¬ 
vidual  layer,  and  imperfect  adhesion  of  these  layers.  Com¬ 
posites  thus  exhibit  additional  modes  of  failure  associated 
with  phase  separations. 


-1- 


Because  of  the  complicated  internal  structure  of  compos¬ 
ite  systems,  both  the  problems  of  stress  analysis  and  failure 
prediction  are  much  more  difficult  than  for  equivalent  single¬ 
phase  specimens.  A  review  of  some  of  the  analyses  which  have 
been  performed  in  connection  with  composite  systems  is  con¬ 
tained  in  Section  II  of  this  report.  One  fact  which  immerges 
very  clearly  from  that  work  is  that  the  stress  field  in  all 
composites  is  truly  three-dimensional  in  character.  Thus, 
even  the  stress  field  in  a  symmetric  laminate  subjected  to  in¬ 
plane  loading,  a  system  to  be  considered  in  this  report,  can¬ 
not  be  accurately  modeled  by  standard  two-dimensional  methods 
of  analysis.  The  review  of  previous  work  in  this  area  further 
indicates  that  relatively  little  work  has  been  carried  out 
with  regard  to  fracture  analysis  of  laminates.  This  fact  is 
clearly  associated  with  the  difficulties  of  extending  fracture 
mechanics  to  composite  structures. 

In  particular,  the  extension  of  fracture  mechanics  to 
laminates  requires  consideration  of  both  three-dimensional 
effects  and  material  interactions.  Hartranft  and  Sih  [1] 
have  developed  an  approximate  three-dimensional  theory  for 
a  single  material  plate  containing  a  through  crack  in  a  re¬ 
mote  tensile  field.  The  approach  to  fracture  of  laminates 
developed  here,  see  Section  III,  is  based  on  an  extension  of 
this  approximate  theory  to  laminates  and  its  application  to 
the  problem  of  a  through  crack  in  a  laminar  plate  subjected 
to  in-plane  loading.  This  crack  represents  a  preliminary 


-2- 


effort  to  model  the  failure  of  composite  plates.  Further 
work  is  clearly  required;  however,  the  relatively  crude  crack 
model  developed  here  already  gives  some  relevant  information 
concerning  the  dependence  of  the  stress  field  near  the  leading 
crack  edge  on  the  geometric  and  material  properties  of  the 
composite  plate. 

The  three-dimensional  character  of  the  stress  field  in 
single  material  plates  containing  cracks,  which  has  been 
studied  by  Hartranft  and  Sih  [1],  is  clearly  relevant  to 
understanding  the  effects  of  flaws  on  the  behavior  of  lami¬ 
nates.  Thus,  an  experimental  study  of  these  three-dimensional 
effects  on  the  ability  of  the  approximate  solution  [1]  to 
model  them  has  been  carried  out  and  is  discussed  in  Section 
IV.  This  work  constitutes  a  necessary  background  to  both 
the  application  of  approximate  theories  of  this  type  to  lami¬ 
nates  and  experimental  programs  for  studying  failure  of  lami¬ 
nates  . 

In  Section  V,  a  self-contained  boundary  value  problem  is 
presented  which  has  particular  application  to  the  delamina¬ 
tion  process  in  layered  composites.  The  effect  of  imperfect 
adhesion  between  layers  of  a  fibrous  composite  is  studied  by 
considering  a  penny-shaped  flaw  in  the  adhesive  layer  con¬ 
necting  two  anisotropic  half-spaces.  The  results  of  the 
analysis  suggest  that  less  load  is  transferred  to  the  adhe¬ 
sive  (thus,  decreasing  the  likelihood  of  delamination)  as  the 

relative  stiffness  of  the  fibrous  layers  is  increased.  This 
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effect  is  amplified  when  the  relative  thickness  of  the  adhe¬ 
sive  layer  is  decreased. 

For  future  work,  the  stress  solutions  for  the  laminated 
plate  presented  in  Section  III  will  be  further  analyzed  and 
used  to  develop  a  criterion  of  fracture.  Preliminary  experi¬ 
ments  have  indicated  that  an  internal  crack  in  a  sandwiched 
layer  of  the  laminate  can  propagate  with  and  without  delami¬ 
nation  depending  upon  the  geometry  and  loading  conditions .  A 
study  on  the  stress  variations  near  the  free  surface  or  inter¬ 
face  will  be  made  since  they  will  play  an  important  role  on 
the  failure  of  the  laminate  composite  weakened  by  flaws. 
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SECTION  II 


MODELING  OP  LAMINATED  PLATES 

Structural  composites  can  be  classified  as  laminated  and 
multiphase.  Multiphase  materials  can  be  further  divided  in¬ 
to  continuous  fiber  reinforced  materials  and  particulate  com¬ 
posites.  Examples  of  fiber  reinforced  composites  are  glass 
fiber  reinforced  plastics  and  boron  fiber  reinforced  plastics 
while  an  example  for  a  particulate  composite  would  be  whisker 
reinforced  metals.  Laminated  composites  consist  of  layers  of 
multiphase  or  homogeneous  materials  bonded  together.  The  gen¬ 
eral  potential  of  laminates  lies  in  the  fact  that  the  plies 
(layers)  of  the  laminate  can  be  oriented  to  meet  both  the 
load  and  direction  requirements  allowing  for  maximum  struc¬ 
tural  efficiency.  An  important  laminated  composite  in  use 
today  is  made  of  layers  of  unidirectional  continuous  fiber 
reinforced  material.  Generally  there  are  two  types  of  these 
laminated  plates:  cross-ply  and  angle-ply.  Cross-ply  lami¬ 
nates  are  those  in  which  the  orthotropic  axes  of  symmetry  in 
each  ply  are  alternately  oriented  at  angles  of  0  and  90  de¬ 
grees  to  the  plate  axes,  while  in  angle-ply  laminates  the  or¬ 
thotropic  axes  of  symmetry  in  each  ply  are  alternately  oriented 
at  arbitrary  angles  of  +0  and  -0  to  the  plate  axes. 

The  analysis  involved  for  composite  systems  can  be  sub¬ 
divided  into  two  major  categories:  analytical  determination 
of  the  mechanical  properties  of  composite  materials,  and  anal¬ 
ysis  of  structures  such  as  laminated  plates.  In  this  section 
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the  theories  predicting  the  elastic  behavior  of  laminated 
plates  are  reviewed  in  detail,  while  a  brief  discussion  of 
the  methods  evaluating  material  constants  is  provided. 

2.1  ELASTIC  PROPERTIES  OP  COMPOSITE  MATERIALS 

Many  investigators  have  been  involved  in  developing 
analytical  methods  by  which  the  lamina  elastic  constants  can 
be  predicted  in  tefrns  of  the  constituent  material  properties 
and  phase  geometry,  i.e.,  the  shape  and  arrangement  of  the 
filaments  and  volume  ratios  of  the  fibers  and  matrix.  For 
example,  the  stiffness  matrix  of  a  particular  lamina  can  be 
considered  as  a  function  of  the  above  parameters,  such  that 


'ij  "  Cij (Ef ,vf jVf ,Em 


,  v  ,  v  ) 

’  mJ  m 


where  E,  v  and  v  are  the  modulus  of  elasticity,  Poisson's 
ratio,  and  volume  fraction  of  the  filament  and  matrix  mate¬ 
rials,  designated  by  the  subscripts  f  and  m,  respectively. 

The  various  theories  which  have  been  proposed  are:  (1) 
netting  analysis,  (2)  strength  of  materials  method,  (3)  self- 
consistent  model,  (4)  variational  method,  (5)  exact  method 
(elasticity),  (6)  statistical  methods,  (7)  finite  element 
technique,  and  (8)  semi-empirical  method.  There  are  several 
basic  assumptions  common  to  these  theories,  namely: 

(1)  The  fibers  are  linearly  elastic  and  generally 
homogeneous . 
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(2)  The  matrix  is  linearly  elastic  and  homogeneous. 

(3)  Fibers  and  matrix  are  free  of  voids. 

(4)  There  is  perfect  bond  at  the  interface  of  the  fibers 
and  matrix. 

(5)  The  ply  is  initially  in  a  stress  free  state. 

(6)  The  fibers  are  regularly  spaced  and  aligned. 

Furthermore,  it  can  be  assumed  that  the  ply  or  lamina  is 
macroscopically  homogeneous,  linearly  elastic  and  generally 
orthotropic  or  transversely  isotropic. 

These  basic  assumptions  lead  to  an  over-simplified  model 
of  the  actual  ply.  The  following  are  some  of  the  reasons: 

(1)  the  filament  arrangement  within  the  ply  is  generally 
random  and  the  filaments  could  be  misaligned,  (2)  a  transi¬ 
tion  region  between  constituents  exists,  e.g.,  lack  of  com¬ 
plete  bond  and  constituents  interaction,  (3)  presence  of 
initial  stresses. 

A  concise  description  and  evaluation  of  these  theories 
along  with  a  comprehensive  list  of  references  is  given  by 
C.  C.  Chamis  and  G.  P.  Sendeckyj  [2].  For  a  unidirectional 
fiber  reinforced  lamina  an  important  point  to  note  is  that 
in  general  the  exact  theories  contribute  a  very  small  cor¬ 
rection  to  calculated  by  the  rule  of  mixture,  i.e.. 


E1  = 


Ervf  +  Vvm 
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where  is  the  modulus  of*  elasticity  of  the  composite  in  the 
direction  of  the  fiber  axis.  However,  moduli  E2  and  G12  are 
rather  sensitive  to  the  geometry  and  constituent  properties. 

2 . 2  COMPOSITE  PLATES 

In  general,  composite  plates  can  be  classified  into  two 
major  categories:  (1)  sandwich  plates,  (2)  anisotropic  homo¬ 
geneous  or  heterogeneous,  single  or  multi-layered  plates. 

Due  to  the  vastness  of  the  field,  a  very  brief  discussion  of 
sandwich  plates  will  be  given,  while  the  attention  will  be 
concentrated  on  the  discussion  of  laminated  plates. 

2.2.1  SANDWICH  PLATES 

These  composite  structures  may  contain  any  number 
of  distinct  homogeneous  layers  separated  by  isotropic  layers. 
The  layers  may  have  different  thicknesses  and  made  of  differ¬ 
ent  materials.  Consider  a  simple  sandwich  plate  consisting 
of  a  core  layer  of  thickness  h  and  two  face  layers  of  thick¬ 
ness  t.  It  is  assumed  that  the  facings  are  relatively  thin, 
i.e.,  t<<h  and  for  most  practical  purposes  they  can  be  treated 
as  membranes,  while  the  core  material  is  considered  to  have 
negligible  flexural  rigidity.  Reissner  [3]  considered  the 
bending  of  such  a  plate,  where,  based  upon  the  above  assump¬ 
tions,  he  formulated  the  problem  such  that  the  transverse 
shears  are  predominantly  carried  by  the  core  plate,  while  the 
bending  stresses  are  primarily  taken  by  the  face  plates.  It 

is  also  assumed  that  the  stresses  in  the  face  plates  are  uni- 
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form  across  the  thickness,  while  no  in-plane  stresses  o  ,  a  , 

*  x5  y5 

t  are  acting  through  the  core  material.  Hence  from  the 

xy 

equations  of  equilibrium  it  can  be  seen  that  the  transverse 
shear  stresses  remain  constant  across  the  thickness  of  the 
core.  Thus  the  strain  energy  of  the  sandwich  plate  can  be 
written  as 


"s  =  §7  !!  [ox,f  +  °y,f  -  2vox,f  •  °y,f 


+  2(1+v)Txy,f]toV  +  jg 


c  -h/2 


x  [t*  „  +  t 2  jdzdxdy 


xz,c  yz  ,c 


(2.2.1. 1) 


where  the  subscripts  f  and  c  refer  to  the  face  layers  and 
core  regions  respectively.  The  relations  between  stresses 
and  couples  are 


ax,f  =  1 


i+tT1  CTy  ,f  =  1  t  (h+t )  *  Txy,f  =  1  t(h+t 


'xz,c  h 


V 

T  =  _Z 
’  yz,c  h 


(2. 2. 1.2) 


Substitution  of  the  above  into  the  strain  energy  relation  re¬ 
sults  in 


ns  =  II 


+  My 


2vM  M  +  2 ( 1+v )M2  ] 
x  y  xy 


+  2G~h  CQx  +  Qy]}dxdy 
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(2. 2. 1.3) 


Note  that  if  G  h  +  ®.  then  the  formulation  is  that  of  the 
c 

12  1 

classical  plate  theory,  where  is  replaced  by  ^ (h+t )'2E 

Hence,  it  can  be  seen  that  the  results  predicted  by  the 
classical  plate  theory  are  more  conservative,  since  it  is 
assumed  that  the  plate  has  infinite  shear  rigidity,  which 
can  result  in  an  unrealistic  solution  when  the  core  layer 
shear  modulus  is  small. 


’  s 


2.2.2  LAMINATED  COMPOSITE  PLATES 


In  this  section,  we  consider  plates  which  are 
made  of  two  or  more  unidirectional  fiber  reinforced  plies 
that  are  bonded  together  to  form  a  laminated  plate.  If  the 
axes  of  orthotropy  of  all  layers  are  oriented  in  the  same 
direction,  then  the  plate  can  be  considered  as  a  unidirection¬ 
al  anisotropic  structure.  If  the  orthotropic  axes  of  symmetry 
in  each  layer  are  oriented  in  a  particular  direction,  then 
such  a  configuration  is  an  angle-ply  laminate.  A  special  case 
of  angle-ply  is  the  cross-ply  laminate,  where  the  orthotropic 
axes  of  symmetry  in  each  ply  are  alternately  oriented  at 
angles  of  0°  and  90°  to  the  plate  axes. 

One  of  the  most  widely  used  theorems  in  formula¬ 
tion  of  various  mathematical  models  of  laminated  plates  is  the 
classical  plate  theory,  where  it  is  assumed,  (1)  the  plane 
sections  normal  to  the  middle  plane  of  the  plate  remain  nor¬ 
mal,  (2)  the  state  of  stress  in  each  layer  is  that  of  gener¬ 
alized  plane  stress  or  plane  stress,  (3)  the  strains  in  the 
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z-direction  are  negligible. 


The  most  elementary  formulation  is  that  of  unidi¬ 
rectional  plane  anisotropic  plate  given  by  Calcote  [4],  where 
it  is  assumed  that  the  small  deformation  theory  is  applicable. 
Using  x,y,z  cartesian  coordinate  system,  and  assuming  linear 
strain  distribution  through  the  thickness 

u  =  -  z  |^,  v  =  -  z  w  =  w(x,y )  (2.2.2. 1) 


then 


ex  =  - 


9  W 

dxT> 


e  = 

y 


3y 


w 

T, 


Y 


xy 


=  -  2z 


9  2w 
9x9y 


(2. 2. 2. 2) 


considering  the  general  case  where  the  principal  material  di¬ 
rections  do  not  correspond  with  the  x  and  y  directions  of  the 
plate,  then  the  stress-strain  relationship  can  be  written  as 


{a}  =  [C]{e} 


(2.2.2. 3) 


where  [C]  is  a  full  3X3  matrix,  because  the  coordinate  planes 
xz  and  yz  are  not  the  planes  of  elastic  symmetry  of  the  lami¬ 
nate.  Defining  stress  couples  in  the  usual  way 


h/2 

Mx  -  _J/2  cVzdz  ■  -  tDllw,xx  +  D12w,yy  +  2Dl6w,xy] 
h/2 

My  =  Jn  °yzdz  ’  -  [D12W,XX  +  D22w,yy  +  2D26w,xy;l 
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where 


h/2 

/ 


-h/2 


x  .  zdz 

xy 


[D,  /"W  +  D0/-w  +  2D/-/-W 

16  , xx  26  ,yy  66  s 


If  the  orthotropic  axes  of  the  material  are  coincident  with 
the  plate  axes,  then  C.^  =  C2g  =  0  implies  =  D 2g  =  0; 
hence  equation  (2.2.2. 4)  reduces  to 
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D,,w  +  2  (D,  0+2D/-/r  )w 

11  ,xxxx  v  12  66'  ,xxyy 


+  Dnow  =  p 

22  ,yyyy  v 


(2. 2. 2. 5) 


which  Is  given  In  [5]  as  the  equation  for  a  plane  orthotropic 

plate.  Further,  If  the  laminate  Is  homogeneous  and  Isotropic 

Eh  ^ 

then  =  (D12+2Dgg)  =  D  =  l-v^ )  *  hence  equation 

(2. 2. 2. 5)  is  reduced  to 


V4w(x,y )  =  g  (2. 2. 2. 6) 

which  is  the  classical  plate  equation. 

Now  let  us  consider  a  laminate  made  of  two  layers 
(angle-ply),  where  two  orthotropic  layers  of  equal  thickness 
are  laminated  in  such  a  way  that  the  axes  of  elastic  symmetry 
form  an  angle  +0  with  the  x,y-axes  in  one  layer  and  an  angle 
-0  in  the  other  layer.  Smith  [6]  formulated  this  problem,  as 
suming  a  linear  strain  distribution  through  the  thickness  as 
given  in  equations  (2. 2. 2.1)  and  (2. 2. 2.2),  and  a  state  of 
generalized  plane  stress.  His  formulation  provided  a  consti¬ 
tutive  relation  exactly  as  equation  (2. 2. 2. 5),  where  he  con¬ 
cluded  that  the  laminate  behaves  as  a  homogeneous  orthotropic 
plate . 


Reissner  and  Stavsky  [7]  pointed  out  this  miscon¬ 
ception  by  reformulating  the  problem,  where  they  considered 

the  effect  of  in-plane  stretching.  This  was  done  by  defining 
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the  strain-displacement  relations  as 


e  =  e  +  zy  ,  e 
x  x  Ax*  y 


e  +  zy  ,  e 

y  Ay  *  xy 


+  ZX, 


(2. 2. 2. 7) 


ex  "  u,x’  ey  v,y’  exy  “  u,y  +  v,x 


(2. 2. 2. 8) 


Xx  w,xx’  xy  w,yy*  Xxy 


(2. 2. 2. 9) 


where  the  mid-surface  strains,  and  bending  curvatures  are 
those  of  the  classical  plate  theory.  The  equilibrium  rela¬ 
tions  are 


N  +  N  +  p  =  0,  N  +  N  +  P  =  0  (2.2.2.10) 

x,x  xy , y  px  *  xy,x  y,y  y 


M  +  M  -  Q  =  0,  M  +  M  -  Q  =  0  (2.2.2.11) 

x,x  xy ,y  ^x  5  xy,x  y,y  ^y 


Q  +  Q  +  P  =  0 
x,x  *y,y 


(2.2.2.12) 


where  the  stress  resultants  are 


(N  ,N  ,N  )  =  (  (a  ,o  , t  )dz 

v  x5  V  xy J  _hJ/2  x5  y5  xy ' 


(2.2.2.13) 


and  the  stress  couples  are 


(Mx’VV  '  J/2  (0x>°y>Txy)zlSz 


(2.2.2.14) 
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The  stress-strain  relationship  is 


{a}  =  [C]{e }  (2.2.2.15) 

Substitution  of  equation  (2.2.2. 7)  into  equation  (2.2.2.15) 
and  the  results  into  equations  (2.2.2.13)  and  (2.2.2.14)  leads 
to 


N 

A 

B 

M 

B 

D 

X 

(2.2.2.16) 


where  A,  B  and  D  are  given  by 


h/2 

( A^ j j j j )  —  J  (IjZjZ  ) • ^ijdz , 

(i,J  =  1,2,6) 


(2.2.2.17) 


Partial  inversion  of  (2.2.2.16)  results  in 


A* 

B* 

N 

M 

C* 

D* 

X 

(2.2.2.18) 


where  A*  =  A-1,  B*  =  A-1.B,  C*  =  B.A"1,  D*  =  D  -  B. A  1.B. 


In  order  to  obtain  the  governing  equations  of 
equilibrium  and  compatibility,  introduce  the  Airy  stress 
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function  <j>(x,y)  such  that  in  the  absence  of  body  forces  the 
stress  resultants  can  be  written  as 


N  =  <J>  ,  N  =  <f> 

V  T  TTTT  ^  TT  “ 


X 


»yy 


XX- 


(2.2.2.19) 


substituting  (2.2.2.19)  and  (2. 2. 2. 9)  into  (2.2.2.18)  and  the 
results  into  the  following  equilibrium  and  compatibility  equa¬ 
tions 


M  +2M  +  M  +  p  =  0 

x  >  xx  xy  ,xy  y,yy  p 


(2.2.2.20) 


°  +  °  _  o 

x,yy  y5xx  “  xy,xy 


(2.2.2.21) 


yields  the  following  system  of  simultaneous  fourth-order  equa¬ 
tions 


LjW  -  L3<j)  =  p  (2.2.2.22) 

L2<J)  +  L^w  =  0  (2.2.2.23) 

where  the  operators  L.  are  given  in  [7]. 

Note  that  equations  (2.2.2.22)  and  (2.2.2.23)  are 
coupled  through  the  operator  L^.  The  coupling  of  w  and  <j> 
enters  into  the  boundary  conditions  as  well,  i.e.,  in  a  prob¬ 
lem  of  transverse  loading,  if  one  has  prescribed  values  of 
M  and  R  =  Q  +  M  at  an  edge  x  =  constant,  at  the  same 

X  XX  XV  •  V 
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time  he  must  formulate  conditions  pertaining  to  N  and  N 
or  u  and  v.  Reissner  and  Stavsky  [7]  also  pointed  out  this 
cross-elasticity  effect  for  stresses,  by  writing  the  stress- 
strain  relationship  in  terms  of  in-plane  strains  and  bending 
curvature  as 


io}  =  [C]{e0}  +  z[C]{X} 


(2.2.2.24) 


Now  invert  equation  (2.2.2.16) 


a  b 


c  d 


N 


M 


(2.2.2.25) 


and  substitute  (2.2.2.25)  into  (2.2.2.24)  to  get 


{a}  =  [C][([a]  +  z[c]){N}  +  ( [b ]  +  z[d]){M}]  (2.2.2.26) 


Thus  it  can  be  seen  that  in  general  a  depends  not  only  on  N 

and  M  but  also  on  N  ,  M  ,  N  ,  M  and  an  analogous  cross 
a  y  y  *y  -x-y 

elasticity  effect  is  encountered  in  the  expressions  for  0 

y 

and  x 

xy 

2.2.3  N- LAYER  LAMINATES 

Dong,  Pister  and  Taylor  [8]  extended  the  above 
formulation  for  multi-layered  anisotropic  plates.  The  as- 


sumptions  of  the  classical  plate  theory  were  used  along  with 
the  strain-displacement  relations  given  in  equation  (2.2.2. 7)* 
Each  layer  and  also  the  laminate  as  a  whole  are  assumed  to  be 
in  a  state  of  generalized  plane  stress.  Calcote  [4]  has  given 
complete  formulation  of  this  problem.  The  stress-strain  re¬ 
lationship  for  the  kth  layer  can  be  written  as 

[ak]  =  [Ck] [e]  (2. 2. 3.1) 


where  [e]  is  given  in  equation  (2. 2. 2. 7),  and  the  symmetric 
matrix  [C^k^]  represents  the  elastic  coefficients  associated 
with  the  kth  lamina.  For  a  particular  layer  these  coeffi¬ 
cients  depend  on  both  the  elastic  properties  of  the  material 
principal  axes  and  the  orientation  of  these  axes  with  respect 
to  the  laminate  axes.  Once  a  reference  surface  for  strains 
and  curvatures  is  defined,  then  stress  resultants  and  stress 
couples  can  be  formulated  with  respect  to  this  surface  as 


h, 
n  k 


(N  ,N  ,N  )  =  l  J  [a;KV  K',x^]dz 
x’  yJ  xy'  l  J  x  y  *  xy 

k  x  I1k_1 


(2. 2. 3. 2) 


n  hk 


.(k)  (k)  ^  (k) 

'x''"y a"xy '  l  5  y 

K  1  nk-l 


«v„)  =  I  /  K  -V  ’Uy  ]zdz  (2 . 2 . 3 •  3) 


Performing  the  above  integration  and  summation. 


N 

A 

_ 

B 

o 

e 

M 

B 

D 

X 

(2. 2. 3. 4) 
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where 


Au  ■  Jj.  cij)(hk  -  hk-i> 

Blj  =  5  Jj  Cij)(hk  "  hk-l> 


D 


ij 


h 


3 

k-1 


) 


Equilibrium  equations  are  the  same  as  for  the  two-ply  laminate 
given  by  equations  (2.2.2.10),  (2.2.2.11)  and  (2.2.2.12),  and 
the  compatibility  equation  (2.2.2.21).  The  governing  differ¬ 
ential  equations  of  the  laminate  can  be  formulated  either  in 
terms  of  an  Airy  stress  function  <f>(x,y)  and  the  transverse  de¬ 
flection  w,  or  in  terms  of  the  three  displacements  u°,  v°  and 
w.  The  system  of  two  simultaneous  fourth-order  equations  ob¬ 
tained  in  terms  of  <j>(x,y)  and  w(x,y)  is  identical  to  that  of 
the  two-ply  laminate  given  in  equations  (2.2.2.22)  and 
(2.2.2.23).  The  formulation  in  terms  of  the  displacements  is 
given  by  Whitney  and  Leissa  [9],  and  Whitney  [10].  The  gov¬ 
erning  equations  for  cross-ply  laminates  are 


L^u°  + 
L2u°  + 

T  O 

-  l4u 


L 2v°  -  L^w  =  0 

L0v°  -  Lr-w  =  0 
3  D 

-  L,-v0  +  L6w  =  p 


(2.2. 3.6) 
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and  the  governing  equations  for  angle-ply  laminates  are 

LjU0  +  L2v°  -  LjW  =  0 

L2u°  +  L3v°  -  LgW  =0  (2.2. 3- 7) 

T  O  t  O  ,  T 

-  L^u  -  LgV  +  LgW  =  p 

where  are  operators  given  in  [9,10].  It  should  be  noted 
that  these  operators  are  different  from  those  given  in  equa¬ 
tions  (2.2.2.22).  P  is  the  transverse  load  acting  on  the 
plate.  These  equations  can  be  simplified  significantly  when 
the  laminated  plate  is  a  symmetric  one  with  respect  to  the  mid¬ 
dle  plane,  i.e.,  when  the  middle  plane  of  the  laminate  is  se¬ 
lected  as  the  reference  surface  of  the  plate  and  the  laminate 
is  composed  of  pairs  of  laminas  such  that  the  two  laminas  of 
any  particular  pair  are  identical  in  thickness  and  orientation. 
Hence,  under  these  conditions  B^.  =  0,  which  implies  no  cou¬ 
pling  between  stretching  and  bending. 

Pagano  [12-14]  has  investigated  the  limitations  of 
the  classical  plate  theory  (CPT)  by  comparing  the  solutions  of 
several  specific  boundary  value  problems  in  this  theory  with 
the  exact  solutions  (within  the  theory  of  linear  elasticity); 
in  particular,  cylindrical  bending  of  laminates  and  transverse 
deflection  of  rectangular  bi-directional  composites. 
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Pagano  [12]  formulated  the  cylindrical  bending 
of  a  laminate  composed  of  m  orthotropic  layers  simply  sup¬ 
ported  on  the  ends  x  =  0  and  x  =  l,  by  assuming  the  plate  to 
be  in  a  state  of  plane  strain.  Since  each  layer  is  ortho- 
tropic,  the  strain-stress  relations  can  be  written  as 

ex  =  Rllcx  +  R13az>  ez  =  R13ax  +  R33°z 
Uz  =  R55txz 

where  are  the  reduced  compliance  coefficients  for  plane 

strain,  defined  in  terms  of  S^.  the  compliances  with  respect 
to  the  principal  axes  of  the  material.  The  equilibrium  equa¬ 
tions  along  with  the  strain-displacement  relations  are  given 
as 


O  +T  =  0  ,  0  +T  =0 

X,X  XZ  ,Z  J  17  17  ’ 


ev  =  u  .  e„  =  w  .  yvr7  =  u  „  +  w 


,XJ  z 


Z  7  '  XZ 


,X 


Note  that  all  the  stresses,  strains  and  displacements  are 
functions  of  x  and  z  only.  A  solution  of  the  following  form 
is  assumed: 


=  f^(z).sinpx 
=  -  p2f .  (z) . sinpx 

✓  •  \ 

=  -  pfj!  (z).COSpX 
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which  satisfies  the  equations  of  equilibrium  identically,  and 
the  eigenvalues  p  =  p(n)  are  determined  from  the  boundary 
conditions.  Substitution  of  these  assumed  solutions  into  the 
stress-strain  and  strain-displacements  relations  results  in  a 
fourth-order  ordinary  differential  equation  in  the  variable 
f^(z)  for  each  layer  which  is  solved  exactly.  In  a  subsequent 
work  Pagano  [13]  solved  the  problem  of  bending  for  a  simply 
supported  rectangular  laminate  under  transverse  load.  Each 
layer  was  assumed  to  be  orthotropic.  A  solution  to  the  gov¬ 
erning  equations  of  elasticity  of  the  form 

u(x,y,z)  =  U(z)cospx. sinqy 

v(x,y,z)  =  V(z ) . sinpx. cosqy 

w(x,y,z)  =  W(z ) sinpx. sinpy 


and 


(U,V,W)  =  (U#,V*,W*).exp(sz) 

was  assumed.  The  eigenvalues  p  and  q  as  well  as  the  constants 
U#,  V*,  and  W*  were  determined  by  substitution  of  the  above 
assumed  solution  into  the  field  equations  and  boundary  con¬ 
ditions  . 


Comparison  of  the  solutions  obtained  by  Pagano 
[12,13]  with  that  of  the  classical  laminate  plate  theory  so¬ 
lutions  discussed  previously  [7-11]  indicate  clearly  that  the 
classical  laminated  plate  theory  (CPT)  underestimates  the 
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plate  deflection  and  gives  a  very  poor  estimate  for  relatively 
low  values  of  S  (span  length/plate  thickness),  e.g.,  In  the 
case  of  cylindrical  bending  at  S  =  20,  the  error  Is  approxi¬ 
mately  20$  and  it  grows  significantly  with  lower  values  of  S; 
however,  for  large  values  of  S  the  solution  approaches  asymp¬ 
totically  the  CPT  results.  Another  important  point  to  note 
is  the  validity  of  the  Kirchhoff  hypothesis,  where  the  solu¬ 
tion  indicates  that  the  assumption  of  linear  displacement 
within  each  layer  may  be  reasonable;  but,  the  deformed  con¬ 
figuration  of  the  original  normal  plane  cannot  be  described 
by  this  hypothesis  for  low  values  of  S.  Meanwhile,  it  can  be 
seen  that  the  deformed  normal  tends  to  straighten  as  the  value 
of  S  increases. 

In  summary,  the  above  comparison  clearly  shows 
the  necessity  of  incorporating  the  influence  of  transverse 
shear  deformation. 

Whitney  [15]  considered  the  bending  of  symmetric 
angle-ply  and  cross-ply  laminates  where  he  included  the  trans¬ 
verse  shear  deformation  in  the  following  way.  Using  the  con¬ 
ventional  x,y,z  coordinate  system,  the  plate  surfaces  z  = 
are  assumed  to  be  free  of  shear  tractions.  Each  layer  is  as¬ 
sumed  to  be  an  orthotropic  material.  Neglecting  the  effect 
of  transverse  normal  stress  a  on  the  gross  response  of  the 
laminate,  the  stress-strain  relation  for  the  kth  layer  can  be 
written  as 
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where  Q  is  the  reduced  stiffness  matrix.  Furthermore,  it 

X  J 

is  assumed  that  the  transverse  shear  stresses  of  the  kth 
layer  have  the  following  form: 


°xz)=  W55)f<z>  +  a55)]'f’x(x>y) 

+  [Qjj|pf(z)  +  ](J>y  (x,y ) 

ayf  =  ^45)f(z)  +  a45)](|)x(x^) 

+  [Q^f(z)  +  aj^]<i>y(x»y) 

In  the  above  ai^  are  constants  which  can  be  determined  from 

X  J 

the  condition  of  the  continuity  of  shear  tractions  at  the  in¬ 
terface  of  adjacent  layers. 

Using  the  strain-displacement  relations  along  with 
the  stress-strain  equations  and  the  assumption  that  w  =  w(x,y) 
yields 
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u(k)  =  -  zw^x  +  [J(z)  +  g1(z)  (k)]<f>x  +  g2(z)  (k)(f.y 
v(k)  =  -  zw ^ y  +  g3(z)(k)(j)x  +  [J(z)  +  g4(z)(k)]<f)y 


where 


J(z)  =  /  f (z )dz 


g^z) 


<k>  -  v(k) 


a<k> 

1 


(i  =  1,2,3, 4) 


and  where  the  constants  d|k^  can  be  determined  from  the  con¬ 
tinuity  condition  of  the  displacements  u  and  v  at  the  inter¬ 
face  of  adjacent  layers,  noting  that  u  =  v  =  0  at  the  mid¬ 
plane  of  the  plate.  These  equations  for  J(z)  and  g^(z)^k^ 
along  with  the  strain-displacement  (linear  elasticity)  rela¬ 
tions  and  Hooke's  law  for  the  kth  layer  will  yield  three 

equations  for  the  stress  couples  M  ,  M  ,  M  as  functions  of 

x  y  xy 

w(x,y),  4>  ( x 3 y )  and  <f>  (x,y)  and  their  derivatives.  Using 

«y 

these  new  equations  along  with  the  conventional  equations  of 
equilibrium  given  in  equations  (2.2.2.11)  and  (2.2.2.12)  will 
result  in  three  constitutive  equations  for  w,  <J>x  and  <j>y. 

These  constitutive  equations  are  rather  formidable;  however, 
for  laminates  in  which  Q^)  =  ajk^  =  0,  they  are  considerably 
simplified  [15].  This  condition  arises  when  the  preferred 
directions  of  the  layers  coincide  with  the  plate  axes. 
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In  order  to  assess  this  formulation,  Whitney  con¬ 
sidered  a  simply  supported  rectangular  laminated  plate  sub¬ 
jected  to  sinusoidal  transverse  loading.  In  order  to  obtain 
relevant  numerical  results,  a  knowledge  of  the  nature  of  the 
function  f(z)  is  necessary.  From  solutions  provided  in  [12] 
and  other  references  mentioned  in  [15] »  Whitney  concluded  that 
the  choice  of 

f(z)  =  1  -  4(|)2 

was  a  judicious  one  for  this  class  of  boundary  value  problems. 
Bending  curves  given  in  [11]  indicate  that  shear  deformation 
increases  the  plate  deflection;  e.g.,  for  a  width  to  depth 
ratio  of  20  this  increase  is  about  15$  for  the  plate  with 
En/E22  =  40,  while  the  laminate  with  E11/E22  =  3  shows  a 
much  less  pronounced  effect,  demonstrating  the  effect  of  an¬ 
isotropy  of  each  individual  layer.  It  can  be  seen  that  de¬ 
spite  the  assumptions  made  on  the  functions  ^xz>  TyZ  an<^  f(z)> 
this  theory  predicts  rather  accurately  the  gross  response  of 
laminated  plates  under  transverse  loads. 

As  is  evident,  bending  theories  of  composite  lami¬ 
nated  plates  have  been  investigated  extensively.  However, 
very  little  attention  has  been  given  to  the  extension  and  in 
particular  the  interlaminar  shear  problems.  Puppo  and  Evensen 
[l6]  have  proposed  the  following  model  for  the  study  of  inter¬ 
laminar  shear  stresses  in  the  laminate.  It  Is  assumed  that 
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the  laminate  can  be  modeled  as  a  set  of  anisotropic  layers 

separated  by  isotropic  layers,  where  the  isotropic  layers 

develop  only  shear  stresses,  acting  as  "adhesives"  between 

the  anisotropic  layers.  It  is  also  assumed  that  the  state  of 

stress  in  each  anisotropic  layer  is  that  of  generalized  plane 

stress,  and  the  in-plane  stresses  a  ,  a  ,  x  and  displace- 

x  y  xy 

ments  u,  v  are  the  averages  of  the  actual  values  through  the 
thickness.  For  a  model  consisting  of  two  layers,  the  equilib¬ 
rium  equations  for  the  anisotropic  layers  are 


<i>  +  Tu>  + ; 


* 

XZ 


X3X 


xy3y 


T(i)  +  0(i)  +  JiZ. 
xy,x  y,y  h. 


=  0 


=  0 


where  the  superscript  (i)  indicates  a  particular  layer,  and 

the  positive  and  negative  signs  correspond  to  i=l  and  i-2  re- 

spectively.  The  interlaminar  shear  stresses  x  and  x  are 

xz  yz 

acting  on  the  faces  of  the  layers,  but  from  the  equilibrium 
point  of  view  it  is  assumed  that  they  are  acting  on  the  middle 
plane  of  each  layer.  This  assumption  is  a  valid  one  if  the 
laminate  is  symmetric  or  warping  is  prevented  by  external  sup¬ 
ports.  The  shear  stresses  in  the  isotropic  layer  are  assumed 
to  be 


Txz  =  °Yxz  =  H(ul  -  u2> 


v  -  “  £(vi  -  v 
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where  u,  v  are  the  displacement  components,  and  G  is  the 

shear  modulus  of  the  isotropic  layer,  while  the  stress-strain 

relation  for  the  anisotropic  material  is  {a}  =  [Q. .]{e)  where 

Q.  .  is  the  reduced  stiffness  matrix.  The  solution  is  obtained 
J 

when  four  functions  u-^,  v-^,  and  axe  found  which  satisfy 
the  constitutive  equations  and  the  particular  boundary  condi¬ 
tions.  For  an  infinite  laminate  under  uniform  unidirectional 
deformation,  it  is  shown  that  the  interlaminar  shear  stresses 
are  zero.  However,  at  the  free  edge  of  a  finite  width  lami¬ 
nate,  the  interlaminar  shear  stress  has  a  finite  value.  Re¬ 
sults  based  on  this  theory  indicate  that  the  magnitude  of  the 

£ 

transverse  shear  x  can  be  approximately  25%  of  the  maximum 

xz 

axial  stress  o  ,  while  in  the  classical  laminate  plate  theory 
these  shear  stresses  are  ignored.  Hence,  it  is  evident  that 
for  a  reasonable  strength  analysis  the  effect  of  the  inter¬ 
laminar  shear  stresses  must  be  included. 

Pipes  and  Pagano  [17]  considered  the  problem  of  a 
finite  width  symmetric  angle-ply  laminate  under  uniform  axial 
strain,  in  order  to  study  the  character  of  the  interlaminar 
shear  stress.  The  stress-strain  relations  for  each  layer  with 
respect  to  the  plate  axes  are 
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and  the  strain-displacement  relations  are  those  of  linear 
elasticity.  The  tractions  are  applied  only  on  the  ends  of 
the  laminate  (x  =  constant)  such  that  the  stress  components 
are  independent  of  x.  Hence  the  equilibrium  equations  are 


'xy  >y 


+  T 


xz  ,z 


=  0,  cr 


y»y 


+  Tyz,z  ~  °s  Tyz,y  +  CTz,z  -  0 


Integrating  strain-displacement  relations  and  enforcing  the 
symmetry  conditions  of  the  displacements  with  respect  to  the 
xy  and  xz  planes  and  the  continuity  conditions  at  the  inter¬ 
faces  yields 


u  =  Cgx  +  U(y  ,z ) ,  v  =  V(y,z),  w  =  W(y,z) 

Since  the  strain  e  is  a  constant,  this  type  of  displacement 

A 

field  is  defined  as  "uniform  axial  extension"  [17].  The  con¬ 
stitutive  equations  of  each  layer  are  obtained  by  combining 
the  stress-strain  relations  with  the  reduced  equilibrium  equa¬ 
tions,  the  strain-displacement  relations  and  the  above  dis¬ 
placement  field.  Prom  the  numerical  results  of  a  four  layer, 
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±45°  symmetric  angle-ply  laminate,  it  can  be  seen  that  the 
in-plane  stresses  converge  to  the  values  predicted  by  the 
classical  plate  theory  away  from  the  plate  boundaries.  How¬ 
ever,  the  interlaminar  shear  stress  t  tends  to  approach  in- 
finity  at  the  free  edge,  indicating  the  possibility  of  the 
presence  of  a  singularity  at  the  intersection  of  the  inter¬ 
face  and  free-edge.  Rybicki  [ 1 8 ]  has  carried  out  a  three- 
dimensional  finite  element  analysis  for  the  problem  considered 
by  Pipes  and  Pagano,  i.e.,  a  laminar  plate  subjected  to  uni¬ 
form  end  displacements.  Rybicki* s  results  for  the  stress 
components  are  in  close  agreement  with  those  of  Pipes  and 
Pagano  except  that  Rybicki *s  interlaminar  shearing  stresses 
do  not  become  singular  at  the  plate  edges.  The  fact  that 
Rybicki  did  not  obtain  a  singularity  for  the  shearing  stresses 
may  be  a  consequence  of  his  finite  element  model  which  does 
not  permit  such  singularities.  The  question  of  whether  these 
shearing  stress  components  attain  a  large  finite  value  or  are 
actually  singular  at  the  plate  edges  is  not  satisfactorily  re¬ 
solved  by  either  the  Pipes  and  Pagano  or  the  Rybicki  calcula¬ 
tion  because  of  the  inherent  approximate  nature  of  both  of 
their  numerical  calculations. 

In  summary  the  classical  laminate  plate  theory 
can  be  used  effectively  when  predicting  gross  behavior  of  thin 
laminates;  however,  it  fails  to  predict  the  interlaminar  shear 
stresses  and  their  effect  on  other  laminate  response.  In  order 
to  predict  the  behavior  of  composite  plates  accurately,  the 
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effect  of  transverse  stresses  must  be  included.  The  mathe¬ 
matical  models  of  Puppo  and  Evensen,  and  Pipes  and  Pagano 
consider  the  effect  of  interlaminar  shear  for  the  special 
case  of  symmetric  laminates.  These  models  can  be  used  as  a 
guide  in  understanding  the  shear  transfer  behavior  between 
layers  of  a  composite  laminate,  but  they  are  inadequate  to 
investigate  the  failure  by  either  delamination  or  fracture 
within  an  individual  layer.  Other  mathematical  models  are 
needed  to  investigate  these  failure  mechanisms  in  conjunction 
with  the  free  edge  effects  due  to  geometrical  discontinuities. 

2.2  A  STRESS  ANALYSIS  FOR  LAMINATES  CONTAINING  FLAWS  - 

A  PROPOSED  THEORY 

The  ultimate  purpose  of  models  and  stress  analy¬ 
ses  for  composites  is  to  enable  strength  analysis  or  the  pre¬ 
diction  of  load  carrying  capacity  for  the  composite  system. 

The  approach  to  be  taken  is  clearly  dependent  on  the  mode  of 
failure  anticipated.  For  laminates,  failure  is  often  asso¬ 
ciated  with  crack  growth  and/or  delamination.  Thus  a  knowl¬ 
edge  of  the  local  effects  of  flaws  or  cracks  on  the  response 
of  laminates  is  needed  to  make  accurate  predictions  of  load 
carrying  capacity.  The  theories  discussed  above  give  a  de¬ 
scription  of  the  stress  field  in  laminar  composites  but  do  not 
take  into  account  the  effects  of  material  imperfections. 

In  the  following  section,  an  approximate,  three- 
dimensional  model  for  laminates  is  developed  and  applied  to 

the  problem  of  a  through  crack  in  a  laminate. 
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SECTION  III 


AN  APPROXIMATE  THREE-DIMENSIONAL  LAMINAR  PLATE  THEORY 

Laminar  composites  have  been  developed  whose  strength  to 
weight  ratio  is  considerably  higher  than  that  of  homogeneous 
materials.  This  had  led  to  a  considerable  interest  in  the 
employment  of  composite  plates  for  structural  application. 
Consequently,  a  great  deal  of  effort  has  been  directed  in  both 
areas,  analysis  and  experiments,  to  determine  the  behavior  of 
such  composite  systems.  In  particular,  it  is  clearly  impor¬ 
tant  to  understand  the  failure  mechanisms  and  develop  crite¬ 
ria  to  predict  failure  conditions.  The  failure  phenomena  for 
composites  is  far  more  complicated  than  for  the  corresponding 
homogeneous  structures  because  fiber-matrix  debonding  and 
interlaminar  separation  often  are  observed  in  conjunction  with 
unstable  crack  growth.  The  present  study  attempts  to  (1)  de¬ 
velop  an  approximate  theory  for  laminar  plates  and  (2)  apply 
this  theory  to  a  laminate  containing  a  through  crack  in  order 
to  predict  the  effect  of  such  a  flaw  on  the  behavior  of  the 
composite  system.  The  analysis  is  particularly  directed  to¬ 
ward  developing  a  better  understanding  of  the  failure  mecha¬ 
nisms,  i.e.,  the  interaction  of  cracks  with  delamination. 

In  the  development  of  approximate  theory  for  laminar 
plates,  it  must  be,  at  the  outset,  realized  that  composite 
plates  contain  a  truly  three-dimensional  distribution  of 
stress  even  under  two-dimensional  loading  conditions.  The 
approach  employed  in  this  work  is  to  extend  the  approximate 

-32- 


three-dimensional  theories  for  homogeneous  plates  developed 
by  Reissner  [3]  and  Hartranft  and  Sih  [1]  to  laminates.  These 
theories  are  based  on  the  variational  principle  of  minimum 
complementary  potential  energy.  An  approximate  form  is  as¬ 
sumed  for  the  stress  field  and  the  complementary  potential 
energy  functional  is  minimized  with  respect  to  admissible 
variations  of  this  field,  i.e.,  those  which  satisfy  (a)  the 
equilibrium  equations  of  three-dimensional  elasticity  and 
(b)  the  traction  boundary  conditions.  The  resulting  solution 
is  approximate  in  that  the  compatibility  or  displacement  con¬ 
tinuity  conditions  are  only  satisfied  in  an  averaged  sense. 
Thus  some  difficulties  can  arise  in  the  application  of  such 
a  theory  to  laminates  because  the  interfacial  conditions  be¬ 
tween  layers  requires  continuity  of  displacements  as  well  as 
tractions.  The  choices  for  the  extension  of  the  three-dimen¬ 
sional  plate  theories  to  laminar  structures  will  be  based, 
in  part,  on  just  this  consideration. 

The  complementary  potential  energy  principle  for  two  (or 
more)  isotropic  material  problems  will  be  developed  and  the 
interfacial  continuity  conditions  will  be  demonstrated.  An 
approximate  stress  field  will  then  be  assumed  and  its  conse¬ 
quences  examined.  This  will  lead  to  an  approximate  theory 
for  laminar  plates  which  will,  in  turn,  be  applied  to  the 
problem  of  a  laminate  containing  a  through  crack  subjected  to 
an  in-plane  loading  situation.  It  should  be  noted  that  the 
same  approach  is  applicable  to  problems  associated  with  bend- 
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ing  of  laminar  plates.  Only  minor  changes  in  the  formulation 
are  needed  to  consider  plate  bending. 

The  approximate  theory  to  be  presented  here  is  capable 
of  dealing  with  laminates  developed  from  homogeneous,  isotrop¬ 
ic  plies.  On  the  other  hand,  more  interest  has  been  expressed 
in  laminates  which  are  manufactured  from  fiber  reinforced 
layers,  thereby  exhibiting  higher  strength  to  weight  ratios. 
Extensions  of  the  present  theory  to  laminates  containing  an¬ 
isotropic  plies,  thus  modeling  the  behavior  of  fibrous  lami¬ 
nates,  are  now  possible  as  a  result  of  the  knowledge  gained 
from  this  idealized  approach. 

3.1  THE  VARIATIONAL  PRINCIPLE  OF  MINIMUM  COMPLEMENTARY 

POTENTIAL  ENERGY  FOR  COMPOSITES 

The  complementary  potential  energy  principle  is  developed 
for  two-material  problems  and  the  results  are  used  to  deter¬ 
mine  the  corresponding  variational  principle  for  laminates. 

Recall,  the  complementary  potential  energy  ($)  of  an  iso¬ 
tropic,  homogeneous,  body  is  defined  as  the  strain  energy  of 
that  body  minus  the  work  done  on  the  portion  of  the  body  sur¬ 
face  (S  )  over  which  displacements  are  specified,  i.e., 

$  =  /  Ydv  -  /  T.u.ds  (3.1*1) 

v  Su 

where  the  strain  energy  density  is 
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Consider  a  body,  made  up  of  two  materials,  which  is  sub¬ 
jected  to  traction  boundary  conditions  along  its  external 
surface.  Figure  3.1.1.  Traction  boundary  conditions  are  as¬ 
sumed  over  the  entire  external  surface  for  convenience  in 
writing  so  that  the  complementary  potential  energy  reduces  to 
the  strain  energy  and  the  only  interfacial  surface  (Sj)  must 
be  considered  separately.  Displacement  boundary  conditions 
can  be  handled  in  the  standard  manner  and  do  not  alter  the 
significant  points  of  the  analysis  to  be  discussed.  The  com¬ 
plementary  potential  energy  of  the  body  described  is  the  sum 
of  the  strain  energies  of  the  portion  made  from  material  (1) 
plus  that  formed  from  material  (2),  i.e.. 


~  1(1 ) 


V)dv 


1(2) 


¥(2)dv 


(3.1.2) 


This  functional  ($)  is  to  be  minimized  with  respect  to  the 
stress  components  subject  to  the  constraints  of 


(1)  equilibrium,  =  0,  p  =  1,2 

1 J  9  J 

(2)  satisfaction  of  traction  boundary  conditions, 
alj)nj  =  T±  on  s£p),  p  =  1,2,  and 


continuity  of  tractions  across  S- 

=  -  oJ?V«> 

!j  J  IJ  j 
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(3) 


/  \ 

where  n.  '  are  the  components  of  the  outward  unit  normal  to 
the  surface  of  the  sub-region  (p).  Thus,  on  Sj, 


The  minimization  of  a  functional  subject  to  a  set  of 
constraints  can  be  accomplished  with  the  use  of  Lagrange 
multipliers.  For  this  case,  the  functional  it  is  introduced 
as 


where  A.^  ,  a|2^,  and  (i  =  1,3)  are  the  Lagrange  multi¬ 
pliers.  The  minimization  of  this  functional  tt  with  respect 
to  admissible  variations  of  the  stress  components  leads  to 
the  set  of  governing  equations  and  material  boundary  condi¬ 
tions.  This  variational  procedure  will  be  carried  out  pri¬ 
marily  to  demonstrate  that  it  leads  to  continuity  conditions 
for  the  displacement  components  across  the  material  interface 
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Sir  =  f  — TTT  dv  +  /  — (4y  6a,^;dv 


V1  9«iJ 


V2  9°iJ 


(2)p/_(2) 


|  6i(nJ1>Solj1)  +  nj2'lSaiJ  )'ds  ”  0 


(3-1.4) 


Application  of  the  divergence  theorem  followed  by  collection 
of  terms  yields 


5"  =  l  (^ry+  x^,6oij)dv 

V1  9<7ij 


+  /  (-rlr  + 


v2  3oij 


(1)  X  o  S„(l)xAl) 


-  J  +  ^)n)±}Sa,, 

o  -L  _L  J 


+  (*{2)  +  61)nj2)6aj2)]ds  =  0 


(3.1.5) 


This  result  can  be  further  simplified  by  recalling  that  on  S- 
+  nj2^cjj)  =  0  which  implies  that 


(2).  (2)  _  _  (1).  (1) 

j  6aij  "  nj  6aij 


Thus 
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=  I  TT7  +  x^j>Soij1)dv 


vi  8oi3 


v2  8oiJ 


-  /  (xj1)  -  xj2))nj1)6a|j ^ds  =  0 


(3.1.6) 


In  order  that  6tt  =  0  for  all  admissible  Sa_jx  ^  in  and  on 

(2) 

and  all  admissible  fiajx  in  Vg ,  it  is  necessary  that 


8oij 


1  +  id)  -  n 

IT  +  h,j  -  0 


in  v, 


^2  (2) 
^Ts7  +  h  j  *  0 

3oij 


111  Vr 


(3-1.7) 


(xj1^  -  xj2))  =0  on  Sj 


Prom  equation  (3-1.1)  it  is  recognized  that  x -  =  and 

8°ij  d 

the  first  two  sets  of  equations  are  the  compatibility  condi¬ 
tions  in  and  v^,  respectively.  They  further  imply  that 


(P)  _ 


-  u^' ,  p  =  1,2,  i  =  1,2,3 


Substitution  of  this  result  into  the  third  set  of  equations 
(3-1.7)  yields 


ud>  =  u'2>  on  S- 
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Thus,  it  has  been  demonstrated  that  the  complementary 
strain  energy  for  a  nonhomogene ous  body  is  the  sum  of  the 
strain  energies  of  the  various  constituents  minus  the  work 
done  over  the  portion  of  the  external  boundary  over  which 
displacements  are  prescribed.  Minimization  of  the  comple¬ 
mentary  potential  energy  with  respect  to  admissible  stress 
fields  leads  to  the  compatibility  conditions  within  each  ma¬ 
terial  and  continuity  of  displacements  across  material  inter 
faces . 

3.2  AN  APPROXIMATE  THREE-DIMENSIONAL  THEORY  FOR  LAMINATES 

The  use  of  variational  principles,  e.g.,  the  principle 
of  minimum  complementary  potential  energy,  to  develop  ap¬ 
proximate  solutions  is  a  standard  technique.  An  approximate 
form  is  assumed  for  the  variables,  in  this  case  the  stress 
components,  and  it  is  substituted  into  the  functional.  Min¬ 
imization  of  that  functional  with  respect  to  admissible  vari 
ations  of  the  assumed  functions  (stress  components)  leads  to 
the  choice  for  those  functions  which  best  approximates  the 
exact  solution  in  an  averaged  sense. 

The  procedure  to  be  employed  for  laminated  plates  is  to 
assume  that  each  of  the  stress  components  in  each  layer  can 
be  approximated  as  a  product  of  a  function  of  the  out-of¬ 
plane  variable  z  multiplied  by  a  function  of  the  in-plane 
variables  x  and  y,  i.e.. 
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no  sum  on  i , j ,p 


(3.2.1) 


0(P>  =  f'?)(Z)g'?>(x,y) 


The  purpose  of  this  assumed  form  for  the  stress  components  is 
to  reduce  the  governing  equations  from  dependence  on  three 
independent  variables  to  dependence  on  two  variables  for 
which  many  standard  solution  techniques  exist.  Substitution 
of  the  assumed  form  for  the  components  into  the  equilibrium 
equations  with  the  additional  condition  that  the  equations 
separate  leads  to  the  more  specific  assumption 


°zp)  =  fp(z)^p)<x»y> 


Ct(p),  t(p)]  =  -  f'U)[Z(p)(x,y),  z^tx.y)] 

xz  3  yz  p  x  30  3  y  30 


(3.2.2) 


[oxP)-  °yP)’  Tx?>;l  *  f;^)[SxP)(x>y)>  SyP)V.y), 
T^p)(x,y)] 


The  conditions  of  traction  continuity  across  material  inter¬ 
faces  are  satisfied  by  taking 


[Z(p)  z(p)  z^P)]  =  [z  ,  Z  ,  Z  ] 
L  x  ’  y  ’  z  J  L  xJ  y  z 


f  (z)  =  f  , (z ) 

p  p+1 


f’(z)  =  f ' . n (z) 

p  p+1' 
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at  the  value  of  z  corresponding  to  the  interface  between 
layers  p  and  p+1. 

It  is  perfectly  reasonable  to  substitute  this  assumed 
form  for  the  stress  field  into  the  complementary  potential 
energy  functional  and  minimize  that  functional  with  respect 
to  f(p)(z),  Zz,  Zx,  Z  ,  S^p),  S^p),  and  T^p)  subject  to  the 
conditions  of  equilibrium  in  order  to  develop  an  approximate 
solution  for  laminar  composites.  It  is,  however,  advantageous 
to  examine  the  resulting  form  of  the  strain  components  so  as 
to  choose  the  approximate  solution  which  comes  closest  to 
satisfying  the  displacement  continuity  conditions  across  the 
material  interfaces.  For  example,  the  strain  component  ex 
is  given  by 

■  I  =  B  [f'(ll(S.-,V  -  vf(z)Zz] 


In  order  for  e  and  thus  u  to  be  continuous  it  is  necessary, 

X  X 

though  not  sufficient,  that  s|p^  and  sfp')  be  the  same  for  all 

x  y 

layers.  Examination  of  the  shearing  strain  component  y 


yields  a  similar  condition  for  T 


(P) 

xy 


xy 

Thus  the  assumed  form 


for  the  stress  field  which  will  be  adopted  is: 


°zP)  =  f(p)(z)Zz(x’y) 

^  3  -  f(p)(z)[Zx(x>y>’  Zy(x-y>]  (3-2-3) 
[«<»,  ajp),  =  f,('p)(z)[Sx(x,y),  Sy(x,y),  Txy(x,y)] 
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where  f  (z)  and  f'(z)  must  match  their  respective  values  for 
P  P 

the  adjacent  layers.  For  notational  purposes  it  is  convenient 
to  rewrite  the  assumed  form  for  the  stress  components  in  an 
equivalent  manner  as : 

a^P)  =  f(z)Zz(x,y) 

[t^,  t^]  =  -  f'(z)[Zx(xJy),  Zy(x,y)] 

CCTxP)j  ayP)j  tx5^  =  f"(z)Csx(xsy) »  Sy(x,y),  Txy(x,y)] 

where  f(z)  must  be  continuous  and  have  a  continuous  first 
derivative  across  the  material  interfaces. 

Substitution  of  this  assumed  form  for  the  stress  compo¬ 
nents,  equation  (3.2.4),  into  the  complementary  potential 
energy  functional  and  the  equilibrium  equations  for  a  laminar 
plate  with  edge  loading  and  traction-free  surfaces  (f  =  f*  =  0 
on  surfaces)  yields 

*  =  !C2E7-{  {f"2[Sj  +  S2  -  2v(p)SxSy 

P  CP)  v(p) 

+  2  (1+v  /  >.  )T2  ]  +  2f  ’ 2  ( 1+v  /  \)(Z2  +  Z2) 

(p)  xy  (p)  x  y' 

-  2v(p)ff"(Zz)(Sx+Sy)  +  f 2Z2 }dv 

-  /  [f"  (Sn^nP)  +  Tns°sP))  "  f^U^lds  (3.2.5) 
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and 


Zx,x 

Sx,x 

Txy,x 


+ 


Z 


y,y 


+  T 

xy,y 


+  s 


y3y 


where 


(3.2.6) 


Sn  -  T  n  =  S  n  +  T  n 
n  x  ns  y  x  x  xy  y 


Snny  +  Tnsnx  =  Txynx  +  Syny 


(3.2.7) 


Z  = 
n 


Z  n  + 
x  x 


Z  n 

y  y 


and  U^p^,  U^p^,  and  are  the  normalj  tangential,  and 

transverse  components  of  the  displacements  prescribed  on  the 
edge  (S^p^)  of  layer  p,  respectively. 

Lagrange  multipliers  (A^)  are  employed  to  insure  satis¬ 
faction  of  the  equilibrium  equations,  and  the  functional 


TT  =  0  -  //  [A-,  (S  +  T  -Z) 
JJ  Ll^x,x  xy,y  x 


+  A„(T  +  S  -  Z  )  +  A_(Z  +  Z 

2  xy ,y  y,y  y'  3  x,x  y,y 


Z  )]dxdy 
(3.2.8) 
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is  formed. 


The  variation  of  it  with  respect  to  Sx,  Syi  TXy>  zx>  Zy3 
Z  ,  and  f(z)  with  the  continuity  conditions  enforced  on  f(z) 

Z 

and  f ' (z )  is 


=  l  l  {f"[2SxSSX  +  2Sy5Sy  -  2v(P)(VSy 

P  V/  N 


+  W  +  4(1+V(p))Txy6Txy]  +  4f,2(1+v(p)  X 


x  (Z  6Z  +  Z  6Z  )  -  2v,  yff"[Z  ( 6S  +  6S  ) 
XX  y  y  (p)  zv  X  y' 


+  (S  +  S„)6Z_]  +  2f 2Z  6Z  }dv  -  /  [f"(6S 

a  y  &  Li  Li  /  \  II  ii 


t"t(p) 


+  6T  -  f’SZ  U^jds] 


-  //  [Xn (6S  +  6T  -6 Z  ) 

J J  L  lv  x,x  xy,y  x' 


+  A0(<$T  +  6S  -  6Z  ) 

2  xy,x  y,y  y' 


+  X0 ( 6Z  +  <SZ  -  6Z  ) ]dxdy 
3  x,x  y,y  z/J  J 


+  Y  —  /  {2[S2  +  S2  -  2v f  sS  S 

2E/  \  ^  L  x  y  (p)  x  y 

P  (P)  V(p)  J  * 

+  2(l+v(p))Txy]f',6f"  +  4(l+v(p)(Z2  +  Zy  )  f '  6f r 


2  v  /  xZ  (S  +  S  )  (f  <Sf "  +  f"6f)  +  2Z2fSf  }dv 
(p)  z  x  .  y  z 


- / 


(p) 


C(s  0 

n  n 


(p) 


+  T  U^)Sf" 
ns  s 


Z  U(p)<5f’]ds] 
n  z 


u 


(3.2.9) 


Note  that  6it  has  been  purposely  separated  Into  two  parts. 

The  first  part  contains  only  variations  of  the  in-plane  func¬ 
tions  S  ,  S  T  ,  Z  .  Z  ,  and  Z  :  while  the  second  set  of 
y  *y  X-  y  2 

terms  includes  only  variations  of  f(z)  and  its  derivatives. 
The  integrals  containing  variations  of  the  in-plane  functions 
can  be  integrated  through  the  thickness  and  those  containing 
variations  of  f(z)  can  be  integrated  in  x  and  y.  Thus  6tt  can 
be  rewritten  as : 


■  I  {I1(SXSSX  +  W  +  2<h  -  WT*y 


A 


+  WSy  +  VV  +  I3<ZXSZX  +  W 


+  I„[Zz(«Sx  +  6Sy)  +  (Sx  +  Sy)«Zz]  +  I5Zz6Zz 


h«Sx,x  +  STxy,y  *  SZx>  -  X2«Txy,x  +  SSy,y  '  «,> 


X3(SZx,x  +  SZy,y  -  5V>dxdy 


-  I  (Un5Sn  +  VTns  +  VZn)d0 

u 


+  l  /P  {a^p'*f"(z)6f"(z)  +  a^f  ’  (z)<5f '  (z) 

^  4- 


p  t 


P-1 


+  a^p)[f(z)6f"(z)  +  f"(z)6f(z)]  +  ajp)f(z)6f(z) 
+  gjp^6f’'(z)  +  $^p)<5f »  (z)}dz 
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(3.2.10) 


where 


I-,  = 


f"2(z)dz 

P  E(P)  t  i 


I.  = 


-  I  tT2^  f  f"2  (z  )dz 

p  e(p)  vi 


2(1+V/  \ )  p 

I,  -  J  F  (p)  /  f ' 2 (z)ds 

3  p  E(p)  Vi 


V  /  \  p 

Ih  =  -  l  /  f(z)f"(z)dz 

*  p  E(p)  Vi 


I  eTT  /  f2(z)az 

P  (p)  tp_1 


un  - 


l  jP  Upp*f"(z)dz 


P  t 


P-1 


u 


~  1  !  U^p)f"(z)dz 


P 


P-1 

t 


u 


=  -  l  I  U(p)f’(z)dz 


P  t 


P-1 


a 


(p)  _ 


=  4y/{  ^  +  Sy  -  2v(p)SxSy  +  2(1+v(p) 


)T 


xy 


a 


(P)  _ 


2 (i+V/  \ ) 

-  //  (Z2  +  Z2 )dxdy 

E(p)  A  X  y 


]dxdy 
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<4P)  -  -  //  MSV  +  S„)dxdy 


(p)  A 


z  X  y  ‘ 


a4p)  -  ett  II  Zztody 

(p)  A 


e.(p)  =  /  (s  u(p)  +  t  u^)dc 

1  inn  ns  s  ' 

Cu 


,<p)  =  -  I  Znu(p)dC 

U 


The  divergence  theorem  and  the  integration-by-parts 
technique  are  applied  to  the  variational  equation  6 tt  =  0 
obtain : 


6*  =  /  r<IlSx  +  I2Sy  +  I4Zz  +  AliX)«Sx 


+  (IlSy  +  Vx  +  +  h,y)SSy 


+  (2(1,  -  I0)T  +  X,  +  A0  )6T 

1  2  xy  l,y  2  ,x  xy 


+  (I3zx  +  \1  +  A3;X)6ZX 


+  (I0Z  +  X0  +  A0  )6Z 

3  y  2  3»y  y 


+  (Ij|(Sx  +  Sy)  +  I5Zz  +  X3)SZz]dxdy 


to 
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-  {  +  VSSn  +  (us  +  V{Tns  +  (uz  +  V*ZnMo 

0 

u 

t 

+  11  [ajp)f,u,(z)  +  (2c^p)  -  a^p))f"(z)  +  ajp)f(z) 

p  t  ,  “■  f 

F  P-1 

+  3^p^  -  g^p^]6f(z)dz  +  £  (a-jp)f"(z)  +  alp^f(z) 

±,zz  2,z  p  1  5 

t 

+  g^)6f’(z)  P  +  l  (-a^f"'(z)  +  otj^f'Cz) 
t  p 

P-1 

t 

-  <4p)f’(z)  -  en(p]  +  3^p))6f(z)  P  (3.2.11) 

Vi 

where 


X  =  X,n,  +  X0n0 
n  11  2  2 

-  -  X,n0  +  X~n0 
S  12  2  2 

and  (n,,^)  are  the  direction  cosines  of  the  normal  to  the 
plate  edge. 

The  governing  equations  and  boundary  conditions  are  now 
apparent,  i.e.,  the  in-plane  variables  can  be  determined  from 
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0 


2(I1  -  VTxy  +  Xl,y  +  X2,x 


I3Zx  +  X1  +  X3,x  =  0 


I3Zy  +  X2  +  X3,y  “  0 


l4(sx  ♦  s  )  +  i5zz  +  x3  -  o 


z  =  z  +  z 

z  x,x  y ,y 


(3*2.12) 


Z  ■  S  +.  T 
x  x ,  x  xy  ,y 


Z 


y 


t  +  s 

xy,x  y3y 


with  Xn  =  -unj  Ag  =  -ug,  and  A^  =  ~uz  on  the  portion  of  the 
boundary  over  which  displacements  are  prescribed.  The  bound¬ 
ary  conditions  associated  with  the  in-plane  variables  S  , 

S  , ...  can  be  either  of  the  traction  type  or  averaged  displace¬ 
ments.  On  the  contour  c  along  the  plate  edge*  one  can  prescribe 
either 


S 

T 

Z 


n 

ns 

n 


or 

or 

or 


X 


n 


X 


s 


-u 


s 


-u 


z 


The  function  f(z) 


for  the  pth  layer  is  governed  by 
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(3.2.13) 


ajp)f""(z)  +  (2aip)  -  a^p))f"(z)  +  ajp)f(z) 


-  o  (P )  fi(p) 

"  P2,z  "  Pl,zz 


with  the  conditions  that,  f(z),  fT(z),  (a|p^f"(z)  +  a^^f 
+  e(p))  and  [a{p)f'"(z)  +  (c^p)  -  a^p))f'(z)  +  6^  -  B^p)] 
must  be  continuous  across  material  Interfaces.  The  traction- 
free  surface  conditions  are  reflected  by  the  requirement  that 
f(z)  and  f'(z)  must  be  zero  on  the  plate  surfaces. 


After  some  nontrivial  manipulations  the  set  of  governing 
differential  equations  can  be  separated  into  the  form: 


Zx  -  V2zx  -  h  (alV2X3  +  a2X3> 

Zy  '  a6v*Zy  =  k  'V!x3  +  a2X3)  (3.2.14) 

a^V4)^  +  a^V2X^  +  a5^3  =  0 


provided  that 


9  Z  9  Z 

93T  +  9^  =  "  a7v2A3  "  a8A3 

r  T4  I2I3  T3  -it 

[(I1+I2)+(I2-I|)  2(I1-I2)]I4  I1 

ai =  [i5d1+i2)+i3ii(-2i5]  +  nprfT 

r  |  I2I3  z3  -|/t  ,j  x 

[(I1+I2)+(I^-I^)  2(I1-I3)](I1+I2) 

a2  =  [I5(I1+I2)+I3I4-2Ijj] 
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=  -  2Iit(I1  -  I2)  +  I1I3 


cl  r-  — 

5 


-  *2> 


'1  2 


a7 


'5 v  1  2 


JjL _ 

•4(I3' 


a8  " 


_JVV 

:5(i1+i2)+i1)d3- 


The  other  in-plane  variables  can  be  expressed  in  terms 
,  Zx,  and  Zy  as 


h  ’  -  Vx  -  <V 


X2  *  -  I3Zy  -  k  <X3> 


1  2 


9x9y 


I,  9Zv  9Z 

_ v  (_JL  +  _JL) 


s  = 

x 


’I  2 


‘1  9x‘ 


■2  9y 


+  [I1I3  -  I4(I1-I2)]  g 


-  Ihh  +  ^ 


V  = 


■1  2 


■i  3y 


'2  9x‘ 


+  [i1i3  -  i4(i1-i2)]  - 


-  tI2I3  +  S^T 


(3.2.15) 


3.3  APPLICATION  OP  THE  APPROXIMATE  THEORY  TO  A  LAMINAR  PLATE 

CONTAINING  A  THROUGH  CRACK 

For  this  analysis,  the  crack  length  is  taken  as  small  in 
comparison  with  the  in-plane  dimensions  of  the  plate.  The 
stress  field  in  the  vicinity  of  the  crack  is  thus  independent 
of  local  boundary  effects  and  the  plate  can  be  modeled  as  in¬ 
finite  in  the  two,  in-plane  directions.  This  assumption  sim¬ 
plifies  the  analysis  significantly  while  maintaining  the  es¬ 
sential  characteristics  of  the  crack  laminar  plate  problem. 

In  particular,  the  problem  of  a  laminar  plate  (symmetric 
about  its  mid-plane)  containing  a  crack  in  a  uniaxial  far- 
field,  stress  state  will  be  considered  (Figure  3*3.1).  The 
crack  is  directed  normal  to  the  loading  direction,  thus  mod¬ 
eling  the  most  damaging  type  of  flaw.  This  problem  can  be 
expressed  as  the  superposition  of  two  auxiliary  problems. 

The  first  auxiliary  problem  is  that  of  an  uncracked  plate 
subjected  to  the  same  loading  conditions  as  the  actual  cracked 
plate.  For  the  second  auxiliary  problem,  the  cracked  plate 
configuration  is  considered  with  no  boundary  loads.  Instead, 
the  crack  faces  are  loaded  with  tractions  of  equal  magnitude 
but  opposite  sense  to  those  found  at  the  crack  location  in 
the  first  problem.  Thus  the  superposition  of  the  two  problems 
yields  the  cracked  plate  with  traction  free  crack  faces  and 
remote  loading  as  originally  described. 

The  first  of  the  auxiliary  problems,  an  uncracked  plate 
in  a  far,  uniform  stress  field,  has  nonsingular  stress  compo- 
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nents  and  as  such  does  not  contribute  to  the  crack-tip  singu¬ 
larity.  Therefore,  the  stress  field  in  the  vicinity  of  the 
crack  tips  for  the  case  of  crack  face  loading  is  identical  to 
that  for  the  remote  loading  situation  of  the  actual  cracked 
plate.  It  is  mathematically  advantageous  to  consider  the 
crack  face  loading  problem  because  it  leads  to  stress  fields 
which  die-out  far  from  the  crack  and  can  thus  be  treated  by 
transform  techniques. 

The  problem  of  a  composite  plate,  made  up  of  isotropic 
layers  distributed  symmetrically  about  its  mid-plane,  con¬ 
taining  a  through  crack  will  be  analyzed  for  the  case  of  sym¬ 
metric  crack  face  loading.  The  loading  is  required  to  be 
symmetric  about  the  plate  mid-plane  as  well  as  the  x  =  0  line 
to  eliminate  any  bending  effects.  The  crack  is  directed  along 
the  x-axis.  Figure  (3-3-1)-  In  general,  the  loading  on  the 
crack  face  will  include  normal  pressure  and  transverse  shear¬ 
ing  tractions.  The  symmetry  conditions  require  that  these 
shearing  tractions  be  zero  on  the  intersection  of  the  mid-plane 
of  the  plate  with  the  crack  faces.  They  also  have  to  die-out 
at  the  top  and  bottom  plate  surfaces.  Thus,  except  in  the 
case  of  very  thick  plates  their  effect  will  be  small  in  com¬ 
parison  with  the  pressure  loading  and  the  crack  face  shearing 
tractions  will  be  omitted  from  this  analysis. 

Consider  the  region  x,y,z  >  0.  The  boundary  and  symmet¬ 
ric  conditions  from  three-dimensional  elasticity  are: 
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on  y  =  0 


x  =  x  =0 
xy  yz 

c  =  p  ( x ,  z ) ,  0<x<a 

y 


(3-3.1) 


u  =  0 

y 


on  x  =  0 


t  —  t  =0 
xy  xz 


u  =  0 

X 


and  on  z  =  0 


T  =  T  =0 

xz  yz 


uz  =  0 


The  equivalent  boundary  conditions  for  the  approximate  theory 


on  y  =  0 


T  =  Z  =  0 
xy  y 

S  =  P(x),  0<x<a 

y  * - 


(3.3.2) 


X2  =  0  ,  x>a 


on  x  =  0 


T  =  Zv  =  0 
xy  x 


x1  =  0 
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and  f(z)  is  required  to  be  symmetric  in  z  to  model  the  sym¬ 
metry  conditions  on  the  z  =  0  plane.  Further  f(z)  and  f’(z) 
are  zero  on  the  plate  surface  z  =  ±h/2. 

The  governing  equations  (3-2.14)  are  to  be  solved  subject 
to  the  boundary  conditions  (3-3-2).  Further,  the  functions 

sx»  Sy,  Txy*  V  Zy*  and  Zz  mus't  die-out  at  large  distances 
from  the  crack. 

3-3-1  SOLUTION  PROCEDURE  FOR  THE  IN-PLANE  VARIATIONS  OF 
THE  STRESS  FIELD 

As  indicated  earlier,  the  solution  procedure  will 
be  based  on  applying  the  Fourier  transform  to  the  x  variable 
to  determine  the  nine  in-plane  variables  S  ,  S  ,...  It  is 

x  y 

convenient  to  recognize  the  symmetry  properties  of  these  vari¬ 
ables  in  advance,  thus  reducing  the  Fourier  transform  to  sine 
and  cosine  transforms  for  the  odd  and  even  functions,  respec¬ 
tively.  It  is  observed  that  the  functions  S  ,  S  ,  Z  ,  Z  , 

x5  y»  y  zs 

X2,  and  X3  are  even  in  x  while  Txy,  Zx,  and  X1  are  odd  in  x. 
Therefore  define 

00 

C  t 

S  ( s,y)  =  /  S  (x,y )cos (sx)dx,  etc. 
x  o  x 

and 

00 

TXy(s,y)  =  /  Txy (x,y )sin(sx)dx,  etc. 

Application  of  Fourier  sine  and  cosine  transforms 
to  the  governing  set  Of  equations  (3-2.4)  leads  to  the  follow- 
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ing  set  of  equations: 


d4  A^ 


3  dy' 


d2ZS 


d2X° 

dy2- 


dy* 


d2Zc 


an  d2Xa 

( s 2  +  i-)Zs  =  s  [ — —  -5-^  +  ( 

v  a6  x  a6  dy 


dy7 


Z  _ 


(S!  +  i-)Z°  = 

a£  y 


ax  d  A^ 


dy 3 


.  4 

'  a4s 

1  & 

|ro 

_  !i 

a6 

a6 

a2 

fl 

a6 

a6 

l5)X3 


s 2 ) X3] 


dXo 
2  >> _ 3 


s2) 


dy 

(3. 3. 1.1) 


provided  that 


dZu 

szs  +  = 

x  dy 


d2A?  c 

a7  dr^  +  (a7s2  “  a8)X3 


The  general  solution  to  equations  (3. 3. 3.1)  which 
satisfies  the  asymptotic  conditions  at  infinity  can  be  written 
in  the  form 


A^Csjy)  =  2Re[A1(s )e-ry] 

Zj(s,y)  =  ^  B(s)e“qy  +  Re[P1(s )e~ry ] 
ZC(s,y)  =  B(s)e_(^y  +  Re[Pp(s)e  ry] 

y 


where 


r  =  ±  [s2 


!i-  ±  i  /!& 

2a3  a3 


)  1 


1/2 


(3. 3. 1.2) 
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q  =  /  s2  + 


[a?-a, (s2-r2 ) ] 

Pl(s)  =  "  2s  ri+igTi^-r^lT'  Al(s) 

[a?-a, (s2-r2 ) ] 

P2(s)  "  "  2r  [l+a6(s2-r2)]  Al(s) 


The  additional  transformed  functions  of  the  in-plane  variables 


c  n  S 


can  be  expressed  in  terms  of  A,,  Z  ,  and  Z  through  equations 

j  x  y 

(3.2.15).  In  particular,  the  functions  Ts  ,  Sc,  and  are 

xy  y  d 

obtained  to  study  the  boundary  conditions: 


T  (s,y)  = 


-1 


Sy(s>y>  -  npcf 


dXQ  I,  * 

+  Zx-sZy>] 


d2X° 

{I1  api+  s2i2x3 


c-i 

+  [I1I3  -  I4(I1-I2)]  ^2L 


(3. 3. 1.3) 


[i2i3  +  i4(i1-i2)]sz|} 


_  p  dA? 

A«(s,y)  =  -  I0ZC  -  -r-± 
2  ,,y  3  y  dy 


The  general  solution  is  now  substituted  into  the 
boundary  conditions  (transformed  appropriately).  The  condi¬ 
tion  Zc(s,o)  =  0  yields 

y 


B(s>  =  Re  (2r  hQXrV]V5>’ 


( 3.3.1 .4) 
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while  the  condition  T  (s,o)  =  0  relates  the  real  and  lmagi- 

xy 

nary  parts  of  A, (s)  as 


Re(2[l  +  ^  <1  -  t^><itag(XrV>:lrAl(5)i  "  0 


It  is  convenient  to  express  this  result  symbolically  as 
Im[A1(s)]  =  $1(s)Re[A1(s)]  where 

0  .  .  .  r  2  ^3_  / (a2“kiai^ (1+a6kl)~ala6k3> 

3x(s)  -  {rR  [s  2ag  (  (i+a6k1)ii  +  (a6k3)2 


+  ri  2i£  k3  C' 


ai+a2a6 


-]}  * 


Xq  al+a2a6 

‘  {rR  2a^  k3  [(l+a6k1)i+(a6k3)i] 


P  2  x3  /(a2-k1a1)(l+a6k1)-a1a6k»i 
ri'-s  2ag  '  ( l+a^k1 ) *+ (a^k3 ) 2 


(3. 3. 1.5) 


kl  ’  2a, 


-  lk3  "  /  (2^ 


_  tSL 
a3 


The  set  of  mixed  boundary  conditions 


S  (x,o)  =  -  p(x) ,  x<a 
«y 


X0(x,o)  =  0  ,  x>a 


remain  to  be  satisfied.  These  conditions  can  be  expressed 
in  terms  of  the  transformed  functions  as 

00 

2  c 

—  /  Sy (s,o)cos(sx)ds  =  -  p(x),  x<a 

(3. 3. 1.6) 

CO 

—  /  Xg(s ,o)cos(sx)ds  =  0,  x>a 
o 

The  functions  Sy(s,o)  and  X2(s,o)  can  be  expressed  in  terms 
of  the  one  unknown  function  Re[A-^(s)]  as 


XgCsjO)  =  f1(s)Re[A1(s) ] 
sJ(s,o)  =  f2(s)Re[A1(s)] 


(3. 3. 1.7) 


where 


f1(s)  =  Re[r(l  +  13)3 

f2(s)  =  -(jT-iyxy  Re{(l+iB)[I1r2  +  I2s2  +  (-qrd^Ig )I3 
+  (r2I1  +  s2I2)I3  +  (s2-r2 )  (I1-I2)Iil]  x 


Thus  the  mixed  boundary  conditions  can  be  expressed  in  the 
form  of  a  set  of  dual  integral  equations  as 


/  g(s)R(s)cos(sx)ds  = 


-  2  p(x) ,  x<a 


/  R(s )cos (sx)ds  =  0,  x>_a 


where 


R(s )  =  f1(s)Re[A1(s)] 


(3- 3-1.9) 


g(s)  = 


f2(s) 


This  set  of  dual  integral  equations  can  be  reduced  to  a 
single  integral  equation  by  the  introduction  of  an  auxiliary 
function,  w(x),  defined  as 


w(x)  = 


—  /  R(s)cos(sx)ds,  x<a 
"  0 


,  x^a 


(3.3.1.10) 


Then  R(s)  =  /  w(x)cos (sx)dx  and  the  set  of  equations  is  satis- 
o 

fied  for  all  functions  w(x)  which  satisfy 


/  {g(s)  /  w(x)cos(sx)dx)cos(sx)ds  =  -  ^  p(x)  (3-3. 1.11) 


This  equation  (3.3.1.11)  contains  the  crack  tip 


singularity  and  corresponding  care  must  be  taken  in  the  solu- 
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tion  procedure.  Exact  three-dimensional  solutions  for  sym¬ 
metric  crack  problems  [19],  indicate  that  on  all  interior 
points  in  the  immediate  vicinity  of  the  leading  crack  edge 
the  local  stress  and  displacement  fields  must  be  of  the  form: 


k-i(z)  a  1  q 

a  =  -  [cos(|0  +  -5-  sin( 6 ) sin (-5-0 ) ] 

X  /2?  2  2  2 

k  ( 2 ) 

a  =  — -  [cos(l-)  -  i  sin(0  )sin(|-0 )  ] 

y  /2?  2  2  2 


az  =  v(ax  +  ay}  (3.3.1.12) 


x  =  sin( 0 ) cos (1-0 ) ] 

xy  /2^  2  2 


T  T 

xz 5  yz 


0(1)  as  p+0 


where  (p,0)  are  local  polar  coordinates  in  the  plane  normal 
to  the  crack  edge.  Figure  (3* 3.1). 

The  auxiliary  function  w(x)  corresponds  to 
X2(x,o)  and,  as  such,  is  expected  to  be  locally  proportional 
to  the  square  root  of  the  distance  from  the  crack  tip.  This 
information  is  used  to  introduce  a  second  auxiliary  function, 
ijj(t),  in  an  attempt  to  remove  the  inherent  singularity  from 
equation  (3- 3.1. 11).  Define  ijj(t)  by  the  integral  equation 

w(x)  =  /  dt,  x<a  (3.3.1.13) 

x  /tz=x"2' 
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Substitution  into  equation  (3-3-l.il)  yields,  after  some 
simplifications , 

a  00  /  v  x 

/ .  t ip ( t )  [/  “i — -  J  (st )sin(sx)ds]dt  =  -  /  p(x)dx,  x<a 
o  o  s  °  o 

(3-3-1-14) 


where  J  is  the  zero  order  Bessel  function  of  the  first  kind, 
o 

If  the  expectations  mentioned  above  are  correct,  equation 

(3-3-1-14)  will  be  a  nonsingular  integral  equation  for  ij>(t) 

* 

which  is  amenable  to  numerical  solution  procedures  .  The 


function  g(s)  can  be  written  in  the  form 


C  +  g*(s) 


(3-3-1-15) 


—  2 

and  C  is  given  in  the  Appendix  where  g*(s)  =  0(s  )  as 

Substitution  into  equation  (3-3-1-14)  and  use  of  the  known 
solution  to  Abel’s  equation  leads  to  the  Fredholm  integral 
equation  (see  reference  [20]  for  details) 


C/tiKt)  +  /  /rpKn)  /trT  /  sg* (s)JQ(sn)J0(st)dsdr) 
o  o 


(3-3-1-16) 


A  posteriori  checks  of  the  numerical  results  do  indeed  indi¬ 
cate  that  ip(t)  is  a  well  behaved  function. 
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This  equation  can  be  reduced  to  standard  form  by  the  follow¬ 
ing  nondimensionalizations 

n  =  a£,  t  =  axj  x  =  aB,  s  =  — 

cl 

and 

$(t)  =  /xifj  (ax) 

Then  equation  ( 3 -  3 • 1 - 16 )  can  be  written  as 
1 

C$(x)  +  /  $(£)KU3x)d£ 

0 

=  -  f  /r  /  dB  (3.3.1.17) 

"  0  /x2-02 

where 


K(5,x)  =  /fT  /  6g*(|)Jo(ec)Jo(0x)de 

The  accuracy  of  the  numerical  solution  procedure  to  be  em¬ 
ployed  is  improved  by  separating  g#(s)  as 

g*(s)  =  s4N2  +  h(s)  (3.3.1.18) 

g 

where  h(s)  =  0(s  )  as  s-*-°°.  The  constants  H  and  N2  are  given 
in  the  Appendix.  The  kernel  K(£,x)  can  be  rewritten  as 
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K(£,x)  =  /fx  {a2HIQ ( £aN)KQ (xaN ) 

00 

+  J  0h(|)Jo(65)Jo(ex)de},  0<?<x 

o  "  ° 

where  I  and  K  are  the  zero-order  modified  Bessel  functions 
o  o 

of  the  first  and  second  kind,  respectively. 

The  case  of  crack  face  loading  which  does  not 
vary  with  x  is  considered,  i.e.,  p(x)  =  p.  The  kernel  is  now 
evaluated  numerically  and  the  integral  equation  is  approxi¬ 
mated  by  a  set  of  algebraic  equations  via  Simpson's  rule. 
Numerical  solutions  for  $(x)  are  thereby  obtained. 

The  original  unknown  function,  R(s),  is  obtained 
by  substitution  into  equations  (3.3-1.13)  and  (3.3.1.10). 
Carrying  out  the  indicated  integrations,  R(s)  can  be  written 
as 


R(s)  =||  p{V(l)J1(as) 

-  /  If  [^]J-,(as£HdU  (3.3.1.19) 

0  ai*  /l 

The  unknown  transformed  functions  Sc,  Sc,  etc.,  can  be  deter- 

x  y 

mined  directly  from  R(s)  by  back  substitution. 

It  is  expected  that  damage  resulting  from  either 
crack  propagation  or  interfacial  delamination  will  originate 
in  the  most  highly  stressed  region,  the  vicinity  of  the  crack 
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tip.  Thus  our  main  objective  is  to  obtain  the  best  approxi¬ 
mate  solution  for  the  near  tip  crack  field.  Referring  to 
equation  (3. 3* 1.19)  it  can  be  shown  that  the  singularity  at 
the  crack  tip  is  associated  with  $(1)  and  that  integral  ex¬ 
pression  in  this  equation  has  no  contribution  to  it.  Because 
an  asymptotic  description  of  the  crack  tip  field  is  sought, 
the  integral  expression  in  equation  (3-3.1.19)  is  dropped  and 

the  functions  ,  S^,  etc.,  are  expressed  in  terms  of  the 
a  y 

simplified  equation  for  R(s).  These  transformed  functions 
are  inverted  asymptotically  for  large  s.  Adopting  polar  co¬ 
ordinates  r,  r1,  r2,  and  0,  01,  e2  as  shown  in  Figure  (3.3.1), 
it  is  found  that 


S  (x,y )*  =  -  — g—  PY(1)C{-2J-  sin[f-(e,+e5)] 

/r1r2  xlr2  ^  1  d 

r 

-  ^  cos [ 0 0  -  ^■(e1+02)]}  +  0(1) 

Sy (x,y)  =  +  — g—  P'K(l)C{i^-  sin[|(e1+e2)] 

*rlr2 

r 

+  “  cos[eQ  -  ^-( 61+e2 ) ] }  +  o(i) 

Txy(x,y)  =  -  — P¥(1)C{?^-  cos[|(01+e2)]}  +  0(1) 
/r, r„  1  2 


1 

e^,  t>13  b2 ,  C,  m2  are  given  in  the  Appendix  A. 


-65- 


e, (k,b0+k0b, ) 

Z  (x ,y )  =  -  — g—  PY(1)  E— — .  3  1  ]  x 

^  2 


x  cos[eQ  -  e1+e2 ) ] }  +  o(i) 


Zx(Xjy)  =  aFT(l)  i-  {[!l(b1k3+b2k1-^)]  x 


x  ~  J^(as )sin(sx)e  syds 


e,  b?  00 

+  2^(kik3+b2kl+a6  ^  yJ;i>s)sln(sx)e 

Zy(x,y)  =  aPt(l)  |i-  (b1k3+b2k1+^|)  x 


x  J  yj  (as )cos (sx)e  s^ds 


(3.3.1.20) 


Note  that  Z  (x,y)  and  Z  (x,y)  are  nonsingular.  Further,  ap- 

x  y 

proaching  the  crack  tip  (a,o),  i.e.,  taking  the  limit  r->-a, 
0-*O ,  the  near  tip  fields  can  be  written  as 

Ko  (z)/a  0n  30n 

ax  =  — - —  [  cos  ( 2~)  -  jr  sin(0^)sin(— 2^)  ]  +  0(1) 

/2r^ 

Ka  (z)/a  0 1  n  30n 

o  =  — .  —  [cos(p^)  +  j  sin(0.,  )sin(— j=-)  ]  +  0(1) 

y  /2r~  ^ 


Ka  (z)/a  1  30b 

— 2—— —  sin(0,  )cos  (— ^-)  ]  +  0(1) 


x^ 


KCm  ( z )  /a  0 , 

a  =  -  [ cos  (p~ )  ]  +  0(1) 

z  /2K[  2 


-66- 


(3-3.1.21) 


T 

xz 


0(1) 


where 


K  =  afdJe^C 
a  (z)  =  Pf"(z) 


and 


aT(z)  =  P(b1k3+b2k1)^f (z) 

Note  that  the  stress  field  in  the  vicinity  of  the 
crack  tip,  as  given  above,  contains  the  usual  square  root 
singularity  in  the  plane  normal  to  the  crack  front.  The  de¬ 
pendence  of  the  stress  components  on  the  in-plane  variable 
(r1,0i)  obtained  here  is  in  complete  agreement  with  the  exact 
three-dimensional  solution  for  a  crack  in  an  elastic  medium 
[19].  Note,  the  results  of  this  approximate  solution  indi¬ 
cate  that  the  transverse  shearing  stress  components  do  not 
exhibit  an  in-plane  singularity,  a  possibility  considered  by 
other  authors  [17]-  The  dependence  of  the  near  tip  stress 
field  on  the  through-the-thickness  coordinate,  z,  will  be 
discussed  in  the  following  section. 
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3.3.2  THROUGH-THE-THICKNESS  VARIATION  OF  THE  NEAR-TIP 
STRESS  FIELD 

Two  possible  approaches  may  be  proposed  for  the 
determination  of  f^z).  The  variational  principle.  Section 
3.2,  yields  equations  and  matching  conditions  for  the  func¬ 
tions  f^z)  which  govern  the  stress  distribution  in  the  thick¬ 
ness  direction  of  the  1th  layer.  Functions  obtained  from 
these  equations  give  the  "best"  choice  for  f±(z)  consistent 
with  the  initial  approximations,  i.e.,  the  function  of  z  which 
on  the  average  over  the  entire  plate  is  closest  to  describing 
the  true  stress  field.  The  possibility  of  using  this  approach 
was  considered  but  an  alternative  direction  was  chosen  in¬ 
stead.  Exact  asymptotic  analysis  [19]  described  earlier  in¬ 
dicates  that  at  all  interior  points  near  the  leading  edge  of 
the  crack,  the  plate  is  in  a  state  of  plane  strain,  i.e., 

o_  =  v(a  +  a  )  (3* 3- 2.1) 

z  x  y 

Within  each  layer,  the  plane  strain  condition  can  be  employed 
to  obtain  f.(z)  in  the  vicinity  of  the  crack  tip.  Supposi¬ 
tion  of  the  requirement  that  the  approximate  solution  should 
satisfy  equation  (3. 3.2.1)  leads  to  the  following  differential 
equations  for  f^(z): 

f£(z)  +  p^fi(z)  =  0  (3. 3. 2. 2) 

where  p^  is  equal  to  the  asymptotic  value  of  the  function 
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-  Z  /v, (S  +S  )  at  the  crack  tip,  i.e., 
z  i  x  y 


(blk3 


+  bg^)] 


It  can  be  shown  in  general  that  the  solutions  to 
equation  (3- 3*2.2)  do  not  possess  sufficient  free  parameters 
to  satisfy  both  of  the  free  surface  conditions  (f  =  f1  =  0). 

Hartranft  and  Sih  [1]  have  used  the  condition  in 
equation  (3. 3*2.1)  to  derive  the  z-distribution  of  the 
stresses  through  a  single  layer  plate  with  a  crack.  Their  re¬ 
sults  in  terms  of  the  function  f(z)  is  Acos(pz);  where  p  is 
the  free  parameter  in  the  theory  and  A  is  arbitrary.  Such  an 
approximate  analytical  solution,  however,  cannot  be  expected 
to  hold  throughout  the  thickness  of  the  plate.  This  is  main¬ 
ly  because  experimental  observation  has  shown  that  there 
exists  a  layer  of  material  near  the  plate  surface  which  be¬ 
haves  very  differently  from  that  in  the  bulk.  Specifically, 
the  fracture  testing  of  through  cracks  in  finite  thickness 
plates  has  evidenced  two  surface  layers  of  the  material  com¬ 
monly  referred  to  as  the  "shear  lips"  that  undergo  a  consid¬ 
erable  amount  of  plastic  deformation.  However,  the  interior 
portion  of  the  fractured  plate  remains  essentially  brittle. 
This  is  analogous  to  the  flow  of  fluids  in  a  pipe  where  the 
viscous  action  dominates  in  a  layer  of  fluid  close  to  the 
pipe  wall  while  the  flow  characteristics  in  the  bulk  can  be 
adequately  described  by  using  the  potential  flow  theory. 
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In  order  to  impose  a  limitation  on  the  analytical 
result  in  [1],  Hartranft  and  Sih  introduced  two  narrow  slices 
of  finite  thickness  (boundary  layers  )  eh/2  close  to  the  plate 
surface,  where  the  variation  e  was  chosen  as  e  =  l/(2+8h/a). 
The  same  approach  will  be  carried  out  here  for  laminated 
plates.  The  non-dimensional  boundary  layer  of  thickness  e 
(to  be  chosen  later)  is  introduced  in  the  outer  most  layers  of 
the  composite  plate.  In  this  thin  strip  the  function  f(z)  is 
chosen  such  that  it  satisfies  the  traction-free  condition  on 
the  exterior  surface  and  equilibrium  in  the  interior  region. 
The  three  layer  symmetric  composite.  Figure  (3.3*1)*  is  chosen 
as  the  simplest  nontrivial  example.  The  numeral  1  is  used  to 
refer  to  the  interior  layer  while  the  outside  layers  are  de¬ 
noted  by  2.  The  solutions  to  equation  (3. 3 .2. 2)  which  satisfy 
continuity  of  f^  and  f j  across  the  material  interfaces  are: 


fx(z) 

f2(z) 


=  Acos(p-^z),  0<_z£h-^ 

=  A{cos  (p1h-L)cos[p2  (z-h-j^)  ] 

P-. 

-  sin(p1h1)sin[p2(z-h1)]}. 


(3- 3*2. 3) 

hl-z-hl+h2 


where 


If - 

There  exists  a  small  but  finite  distance  to  the  boundary 
within  which  the  resolution  of  experimental  measurements  and 
accuracy  of  any  numerical  solution  breaks  down. 
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(3. 3. 2. 4) 


(b1k3  +  b2k1)] 


P 


2 

2 


P 


2 

1 


See  Appendix  A  for  b-^,  b2,  e^»  mg,  k-^  and  k^. 

The  functions  f^(z)  given  by  equation  (3. 3* 2. 3) 
which  satisfy  the  plane  strain  conditions  are  not  general 
enough  to  allow  imposition  of  the  free  surface  conditions  at 
z  =  h-^  +  h2-  It  is  reasonable  to  introduce  a  boundary  layer 
adjacent  to  the  free  surface  within  which  an  arbitrary  func¬ 
tion  f2(z)  is  constructed  such  that  the  traction-free  condi¬ 
tions  f2(h^+h2)  =  f^(h^+h2)  =  0  are  satisfied  along  with  the 
requirements  that  f2(z)  and  its  first  and  second  derivatives 
must  be  continuous  at  z  =  h^  +  h2(l-e).  This  arbitrary  func¬ 
tion  can  be  represented  by  the  following  fourth  order  poly¬ 
nomial 


f2(z)  =  a1(z-h-^-h2)4  +  a2(z-h^-h2)3  +  a3(z-h^-h2)2 

+  a^(z-h-^-h2)  +  a^,  (3- 3-2. 5) 

hx  +  h2(l-e)  <  z  <  h1+h2 

The  five  coefficients  (a^,  i  =  1,2,.., 5)  are  determined  from 
the  conditions  mentioned  before. 

At  this  point  the  functions  f^(z)  (i  =  1,2)  can 
be  expressed  solely  in  terms  of  the  parameter  p^,  which  de- 
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pends  on  the  behavior  of  the  functions  S  ,  S  ,  and  Z  (x,y). 

The  dependency  of  f^(z)  on  p^  is  found  by  using  an  iterative 

procedure,  which  consists  of  selecting  an  initial  value  of  p^ 

that  specifies  a  given  set  of  the  coefficients  of  asymptotic 

behavior  of  S„,  S,.  and  Z  _(x,y).  This  leads  to  a  new  value  of 
x  y  z 

p^  obtained  from  equation  (3-3-2.1)  and  the  procedure  is  re¬ 
peated  until  it  converges  to  any  one  of  a  set  of  discrete 
values  for  p.j^  which  satisfy  equation  (3. 3 .2.2).  This  itera¬ 
tive  procedure  should  be  modified  such  that  the  solution  for 
p^  converges  to  that  value  which  best  approximates  the  actual 
load  applied  on  the  crack  surface. 

The  crack  tip  stress  field  as  given  by  equations 
(3.3*1.21)  is  now  completely  determined  and  there  remains 
simply  the  numerical  calculation  of  the  stress  intensity  fac¬ 
tor,  k^(z)  which  can  be  expressed  in  the  form 

kx(z)  =  KCp (z )/a 

The  newly  defined  function 

op(z)  =  PfJ(z),  (i  =  1,2) 

represents  the  equivalent  load  distribution  on  the  crack 
faces  which  is  compatible  with  the  variational  theory.  To 
recapitulate,  the  in-plane  stress  components  are  proportional 
to  f£(z),  and  the  transverse  normal  stress  a z  is  equal  to 
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v1(ax+0y)  except  in  the  boundary  layer,  while  the  transverse 
shear  stresses  x  and  x  are  non-singular.  Thus,  the  singu- 

aZ  y  ci 

lar  stress  field  in  the  vicinity  of  the  crack  tip  can  be 
characterized  by  the  non-dimensional  stress  intensity  param¬ 
eter  K  and  the  eigenvalues  p1  and  p2  which  govern  the  func¬ 
tions  f^z). 

It  should  be  pointed  out  that  the  behavior  of 
the  laminate  may  be  qualitatively  different  in  the  vicinity 
of  the  material  interfaces  than  that  within  each  layer.  Such 
an  effect  may  have  significant  influence  on  the  load  transfer 
from  one  layer  to  the  next.  This  line  of  reasoning  suggests 
the  possibility  of  introducing  additional  boundary  layers  or 
transition  regions  at  the  material  interfaces.  However,  the 
introduction  of  such  boundary  layers  is  not  required  by  this 
formulation,  because  the  functions  f1(z)  given  in  equation 
(3- 3 *2. 3)  do  satisfy  the  traction  continuity  conditions  at 
the  interfaces. 

In  order  to  carry  out  the  numerical  analysis  it 
is  necessary  to  assume  a  value  for  e,  the  non-dimensional 
boundary  layer  thickness.  For  laminates  the  boundary  layer 
thickness  will  in  general  depend  on  the  elastic  constants 
Ei/E2,  v2  and  geometric  parameters  h2/h1  and  (h^+h2)/a, 

but  sufficient  data  for  modeling  the  functional  dependence  is 
not  yet  available.  Hence  as  a  first  attempt  in  modeling  lami¬ 
nates  the  present  formulation  includes  boundary  layers  at  the 
exterior  surfaces  but  no  interfacial  transition  layers.  This 
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non-dimensional  boundary  layer  thickness  was  assumed  to  be 
of  the  form 


e  = 


h-,  +h0 

2+16 


hp 


(3. 3. 2. 6) 


The  consequence  of  this  choice  will  be  discussed  in  the  fol¬ 
lowing  section. 


3.3.3  NUMERICAL  RESULTS  AND  DISCUSSION 

Numerical  results  are  obtained  for  the  dependence 
of  the  crack-tip  stress  field  on  the  geometric  and  material 
parameters  (h^/h-^,  h-^/a,  E-^/E^,  v1  and  Vg)  of  the  layered 
plate  containing  a  through  crack  for  two  formulations.  The 
first  does  not  and  the  second  include  the  boundary  layers 
discussed  in  Section  3- 3*2.  The  amplitude  of  the  singular 
stress  field,  equation  (3.3.1.21),  is  presented  in  the  non- 
dimensional  form 


K  =  k1(z)/(op(z)/a) 

where  0  (z)  =  PfV(z)  is  the  crack  face  load  distribution.  The 
P 

results  for  K  as  a  function  of  the  composite  plate  parameters 
are  presented  graphically. 

Results  for  the  formulation  without  any  boundary 
layers  are  presented  first  in  Figures  3-3. 3-1  -  3 .3. 3. 4.  In 
this  case  no  attempt  has  been  made  to  satisfy  the  boundary 
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conditions  on  the  plate  surfaces  z  =  ±(h-^+h2).  As  a  conse¬ 
quence,  the  parameter  p^  which  governs  the  z  dependence  of 
the  stress  field  given  by  equation  C3-3-2.3)  can  be  chosen 
arbitrarily.  Note,  however,  that  for  p^(h1+hg)  greater  than 
tt/2,  the  crack  face  load  distribution  changes  from  compression 


TT 


to  tension  at  z  =  7^—  for  =  v 2  and  no  longer  adequately 
models  the  physical  problem  under  consideration.  For  this 
reason,  the  results  are  confined  to  the  range  0<p1<ir/2. 


Referring  to  Figures  3- 3. 3-1  -  3- 3. 3-4,  and  de¬ 
fining  the  symbol  =  p-^h-^,  note  that  the  stress  intensifi¬ 
cation  parameter  K  is  strongly  influenced  by  the  value  of  p-^ 
chosen,  i.e.,  as  p^  is  decreased,  K  increases  for  finite  h^/a. 
In  the  limit  h-^/a  becomes  large,  K  approaches  unity  for  all 
values  of  p^.  The  reason  for  the  p^  dependence  is  clear  when 
one  recalls  the  z  dependence  of  the  near  tip  field,  equation 
(3. 3-2. 3)*  As  p-^  decreases  the  rate  of  change  of  this  field 
with  respect  to  z  is  correspondingly  decreased.  In  fact,  as 
p^  approaches  zero,  the  crack  tip  stress  field  becomes  inde¬ 
pendent  of  z  and  reduces  to  the  two-dimensional  plane  strain 
solution  for  which  K  equals  1.0. 


The  dependence  of  K  on  the  relative  layer  thick¬ 
ness  h2/h1  is  also  seen  in  these  graphs.  As  h^/h^  increases, 
the  stress  intensification  parameter,  K,  correspondingly  in¬ 
creases  for  finite  h^/a  approaching  1.0  for  large  h2/h^. 

This  result  is  directly  associated  with  the  increase  in  plate 


-75- 


thickness  and  the  approach  of  the  plane  strain  condition  at 
large  h2/h1.  The  dependence  of  K  on  h2/h1  is  more  pronounced 
at  low  values  of  p^. 

Figures  3.3. 3.1  -  3 • 3 - 3 - ^  show  that  the  depen¬ 
dence  of  K  on  the  relative  layer  stiffness  E^/E 2  is  rather 
weak.  Results  for  the  formulation  which  includes  boundary 
layers  and  satisfies  the  mathematical  boundary  conditions  are 
to  be  presented  next.  They  exhibit  a  stronger  K  dependence 
on  E]/E2. 

Boundary  layers  at  the  free  surfaces  are  included 
in  the  refined  formulation  in  order  to  (a)  satisfy  the  mathe¬ 
matical  boundary  conditions  and  (b)  model  the  physically  ob¬ 
served  surface  phenomena  referred  to  as  the  shear  lips  zone 
as  discussed  in  Section  3-3.2.  In  this  formulation  the  pa¬ 
rameter  px  is  no  longer  arbitrary.  Numerical  results  for 

V1  =  v2  (for  which  p2  =  pl^  suSSest  that  pi^hi+h2^  ls  almost 
constant  and  approximately  equal  to  tt/2.  Deviations  are  at¬ 
tributed  to  the  specific  form  of  the  boundary  layer  solution. 
Thus ,  except  near  the  free  surfaces 


f(z) 


/TT  Z 

C0S  2  h1+h2 


) 


This  implies  that  near  the  crack  tip  the  normal  stress  compo¬ 
nents  0,0,  and  0  attain  a  maximum  at  the  center  plane  of 
x  y  z 

the  composite  plate  and  approach  zero  near  the  top  and  bottom 
surfaces  of  the  plate. 
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Numerical  results  are  obtained  for  the  dependence 
of  the  stress  intensification  factor  K  on  the  geometric  and 
material  properties  of  the  laminate  associated  with  the  re¬ 
fined  formulation  which  includes  boundary  layers.  Results 
are  presented  graphically  in  the  form  of  plots  for  K  as  a 
function  of  h^/a  at  various  values  of  E^/E^,  h^/h^  with 
V1  =  v2  =  .3  in  Figures  3  -  3  -  3  *  5  -  3* 3.3-14. 

Attention  is  directed  to  the  influence  of  the 
relative  stiffness  of  the  inside  layer  to  the  outer  layers 
(E^/Eg)  on  the  amplitude  of  the  stress  field  in  the  vicinity 
of  the  crack  tip.  Figures  3* 3. 3-5  -  3-3. 3- 9  indicate  that  K 
decreases  as  E-^/E^  increases  if  the  other  parameters  are  held 
constant.  As  E^/E2  increases  the  load  is  carried  increasingly 
by  the  middle  layer.  In  the  limit  as  E^  goes  to  zero,  the 
actual  plate  half  thickness  is  reduced  from  h^+h2  to  hi  and 

the  stress  intensity  factor  is  expected  to  decrease  according- 

* 

ly  .  The  influence  of  h^/h^  on  this  effect  can  be  seen  by 
comparing  Figures  3 - 3 -  3 - 5  -  3 • 3 • 3  -  9 *  As  h2/h^  increases,  the 
effect  of  changes  of  E^/E2  on  K  becomes  less  pronounced.  To 
understand  these  interactions,  recall  that  at  finite  E-^/E^ 

¥ - 

In  the  present  model,  it  is  not  possible  to  account  com¬ 
pletely  for  this  effect.  The  load  distribution  permitted  in 
this  formulation  for  finite  but  small  E 2  is  qualitatively  dif¬ 
ferent  than  that  for  E 2  =  0,  i.e.,  for  E 2  finite, 

P1  *  2  ^h1+h2^  whlle  at  E2  =  0,  px  s  ^ 
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increasing  h2/h-^  causes  the  outer  layers  to  carry  a  larger 
portion  of  the  applied  load  and  for  very  large  h2/h-^  the 
laminate  behavior  should  be  almost  independent  of  the  mate¬ 
rial  properties  of  the  inside  layer. 

Figures  3* 3 *3*10  -  3* 3* 3*13  explicitly  demonstrate 
the  influence  of  the  parameter  h^/h^  on  K.  If  h-^/a  is  held 
constant  increasing  hg/h-^  increases  the  total  plate  thickness 
and  K  is  observed  to  increase  accordingly. 

Figure  3 -3. 3*14  is  a  plot  of  K  against  h-^/a  for 

the  special  condition  in  which  the  total  plate  thickness  to 

h-,+h2 

crack  length  is  held  constant  ( — - —  =  3.0).  Note  that  for 

a. 

finite  l-^/h^,  K  decreases  as  E^/E^  increases.  In  this  case, 
load  is  transferred  to  the  middle  layer  and  the  "apparent" 
plate  thickness  is  decreased.  Calculations  were  also  carried 
out  for  the  special  case  E^/Eg  =  1.0.  For  E^/E2  =  1.0,  the 
stress  intensification  parameter  K  is  not  expected  to  vary 
with  h-L/h2  as  long  as  the  total  plate  thickness  to  crack 
length  ratio  is  held  constant.  The  slight  variation  observed 
is  associated  with  the  boundary  layer  thickness  dependence  on 
h^/h2  (equation  (3.3 .2.6))  which  is  inappropriate  for  this 
special  case. 

In  summary,  the  approximate  three-dimensional 
stress  analysis  presented  here  for  a  laminate  containing  a 
through  crack  indicates  that  the  amplitude  of  the  near  tip 
stress  field  tends  to  increase  with  plate  thickness.  The 
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stress  intensification  parameter  K  was  found  to  decrease  as 
the  relative  stiffness  of  interior  layer  to  outer  layers 

is  increased.  Poisson’s  ratio  influences  are  to  be 
studied  in  subsequent  work. 

3.3.4  LIMITATIONS  OF  PRESENT  THEORY  AND  FUTURE  WORK 

The  present  formulation  includes  boundary  layers 
at  the  top  and  bottom  surfaces  of  the  composite  plate;  how¬ 
ever,  it  does  not  allow  for  a  transition  layer  at  the  material 
interfaces.  Such  an  interfacial  boundary  layer  is  not  re¬ 
quired  by  this  variational  formulation  as  all  equilibrium 
and  traction  conditions  can  be  satisfied  without  the  intro¬ 
duction  of  interfacial  layers.  However,  the  effect  of  not 
including  these  transition  layers  in  this  approximate  formu¬ 
lation  is  that  the  crack  face  load  distribution  is  fixed  in 
the  z  direction  independent  of  the  material  properties  of  the 
layers.  The  physics  of  the  problem  suggests  that  the  major 
portion  of  the  load  per  unit  thickness  should  be  transmitted 
in  the  stiffer  layers.  Interfacial  transition  layers  are 
needed  in  the  approximate  theory  to  take  into  account  this 
type  of  load  distribution  -  layer  property  interaction.  Ex¬ 
tensions  of  the  present  model  to  include  these  layers  are 
currently  under  study.  A  stronger  dependence  of  the  stress 
distribution  in  the  laminate  on  the  material  properties  of 
the  constituents  is  expected  in  the  improved  model. 

Further  theoretical  and  experimental  work  is  also 

needed  to  employ  the  approximate  stress  analyses  as  presented 
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here  for  the  prediction  of  laminate  strength.  The  possibility 
of  defining  equivalent  fracture  toughness  as  a  failure  crite¬ 
rion  for  composite  plates  is  being  considered.  Some  related 
work  for  finite  thickness  homogeneous  plates  containing 
through  cracks  is  discussed  in  Section  IV  of  this  report. 
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SECTION  IV 


EXPERIMENTAL  MEASUREMENT  OF  THE  EFFECT  OF  PLATE 
THICKNESS  ON  THE  STRAIN  ENERGY  RELEASE  RATE 

The  results  of  an  experimental  program  to  determine  the 
effects  of  plate  thickness  on  the  failure  of  homogeneous 
plates  containing  through  cracks  subjected  to  tensile  loading 
are  presented  here.  This  work  is  a  prerequisite  to  experi¬ 
mental  determination  of  strength  for  laminated  plates  con¬ 
taining  flaws  where  a  large  number  of  material  and  geometric 
parameters  will  have  to  be  considered. 

In  this  study,  compliance  measurements  and  critical  load 
tests  were  performed  for  a  set  of  specimens  covering  a  range 
of  plate  thicknesses  as  well  as  crack  lengths.  The  results 
of  these  tests  enable  an  indirect  calculation  for  the  average 
through- the-thickness  of  the  energy  release  rate;  i.e.,  the 
rate  at  which  the  strain  energy  of  the  plate  is  decreased  as 
the  crack  extends. 

The  energy  release  rate  (or  the  stress  intensity  factor 
which  is  directly  related  to  it)  is  the  basic  parameter  used 

for  failure  prediction  in  two-dimensional  elastic  fracture 

* 

mechanics  .  However,  for  finite  thickness  specimens,  the 
averaged  energy  release  rate  is  expected  to  depend  on  the 

1  - - 

Reference  is  made  to  the  plane  stress  and  plane  strain  formu¬ 
lations  of  two-dimensional  elasticity. 
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ratio  of  the  plate  thickness  to  the  crack  length.  An  approx¬ 
imate  stress  analysis  which  takes  this  thickness  effect  into 
account  has  been  carried  out  by  Hartranft  and  Sih  [1,28], 
(Figure  4.1)  and  is  extended  to  laminates  in  Section  III  of 
this  report.  The  experimental  results  obtained  here  will  be 
compared  to  the  predictions  of  that  approximate,  three-dimen¬ 
sional  model.  An  analogous  theory  has  been  developed  [28] 
and  a  corresponding  program  carried  out  [29]  for  the  bending 
of  finite  thickness  plates  with  through  cracks. 

4.1  CALCULATION  OF  ENERGY  RELEASE  RATE  FROM  EXPERIMENTAL 

RESULTS 

Let  U  be  the  strain  energy  of  the  cracked  plate.  If  the 
crack  extends  under  constant  displacement  conditions  so  that 
no  external  work  is  done  to  the  plate,  then  the  energy  re¬ 
lease  rate  per  unit  of  new  crack  surface  area  A  created  is 
defined  as 

G  =-■—  (4.1.1) 

For  the  present  experimental  investigation,  the  test 
specimen  is  loaded  in  tension  by  concentrated  loads  P  as 
shown  in  Figure  4.2.  The  points  of  application  of  the  loads 
are  displaced  a  distance  8  =  6( P)  during  loading.  Thus  the 
strain  energy  stored  in  the  specimen  as  a  result  of  the  loads 
P  is 
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(4.1.2) 


provided  the  plate  response  is  completely  elastic.  If  the 
crack  in  the  loaded  plate  then  extends  under  fixed  grip  con¬ 
ditions,  i.e.,  if  5  is  held  constant  during  crack  growth, 
then 


G 


9U 

9A 


constant 


(4.1.3) 


At  this  point,  it  is  convenient  to  define  an  apparent 
plate  stiffness  k  =  P/6.  Then 


96  _  p  9 (1/k)  .  1  9P 
9A  r  9A  k  9A 


Thus 


9P 

9A 


constant 


kp  mau 


Substitution  of  this  result  into  equation  (4.1.3)  gives 


G 


1  p 2  9(6/P) 

2  9A 


(4.1.4) 


Thus,  the  critical  energy  release  rate  (the  value  at  fracture 
initiation)  can  be  determined  from  experimental  results  for 
the  load  P  at  fracture  and  the  specimen  compliance  6/P  as  a 
function  of  crack  surface  area. 
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4.2  THE  EXPERIMENTAL  PROGRAM 


Plexiglas  was  chosen  as  the  material  for  the  test  speci¬ 
mens.  It  possesses  the  following  desirable  characteristics: 

a.  Plexiglas  exhibits  comparatively  little  plasticity 
up  to  failure. 

b.  A  cost  material,  plexiglas  is  very  nearly  isotropic. 

c.  Since  plexiglas  is  transparent,  the  detection  of 
crack  extension  in  it  is  relatively  easy. 

The  specimens  were  16  inches  long  and  3  inches  wide. 
Their  thickness  varied  from  1/16  to  3/4  inch.  After  the 
specimens  were  cut  from  a  sheet  of  cast  material,  two  holes 
were  drilled  in  order  to  mount  them  on  the  loading  pins. 
Another  small  hole  of  approximately  1/8  inch  in  diameter  was 
drilled  in  the  center  of  the  specimen  so  that  a  jeweler’s  saw 
could  be  used  to  cut  the  crack.  Twelve  specimens  for  each 
thickness  were  prepared,  (see  Figure  4.3)* 

4.3  TESTING  PROCEDURE 

An  Instron  universal  testing  machine  with  10,000  lb 
capacity  was  used  to  apply  the  tensile  load  between  the  grips. 
A  specimen  was  tested  with  the  small  hole  of  1/8  inch  radius 
without  a  crack  to  determine  the  load-displacement  relation¬ 
ship  for  zero  crack  length.  The  displacements  in  each  test 
were  measured  between  the  points  of  load  application  and  also 
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for  an  8  inch  gage  length,  4  inches  from  each  side  of  the 
crack.  Those  two  different  methods  of  measuring  the  displace¬ 
ment  were  used  in  order  to  study  the  importance  of  the  dis¬ 
placement  gage  length.  Each  specimen  was  accurately  -centered 
before  loading.  The  load  was  applied  in  small  increments. 

At  each  increment  the  relative  displacement  between  the  two 
points  in  consideration  (points  of  load  application  or  8  inch 
gage  length)  was  measured  with  2  dial  indicators  read  to  an 
accuracy  of  ±0.001  inches.  During  the  load  application,  the 
displacement  had  a  constant  rate  of  0.02  inches  per  minute. 

For  each  specimen,  the  loading  was  applied  five  or  more  times 
until  a  repeatable  load-displacement  curve  was  obtained. 

After  determining  the  load-displacement  curve  for  the  speci¬ 
men  with  the  small  hole,  a  small  crack  of  1/4  inch  length  was 
cut  with  a  jeweler’s  saw.  The  specimen  was  then  mounted  on 
the  machine  and  subjected  to  tension,  applying  the  load  in 
small  increments  and  measuring  the  displacements  for  each 
load  increment. 

Following  this  testing  procedure  in  which  the  load-dis¬ 
placement  curve  was  determined,  the  crack  was  lengthened  by 
an  increment  of  1/4  inch  and  tested  again.  The  procedure  was 
continued  until  the  crack  length  of  1  1/2  inches  was  reached. 
This  is  the  maximum  crack  length  advisable  for  avoiding  ef¬ 
fect  of  finite  width  [30].  Details  of  the  test  apparatus 
are  shown  in  Figures  4.4  and  4.5.  The  same  procedure  was  per¬ 
formed  for  each  different  thickness.  Figure  4.6  shows  typical 

load-displacement  curves.  A  second  series  of  tests  were  car- 
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ried  out  to  determine  the  critical  loads  for  each  plate 
thickness  and  crack  length.  In  Figure  4.7,  the  critical  load 
is  plotted  against  crack  length  for  the  case  of  3/8  inch 
thick  specimens . 

4.4  EXPERIMENTALLY  DETERMINED  ENERGY  RELEASE  RATES 

Table  4.1  contains  a  listing  of  the  experimentally  deter¬ 
mined  compliances  (6/P)  and  fracture  loads  (P  )  for  the 
specimens  tested.  For  each  thickness  considered,  compliance 
was  plotted  against  crack  length  and  corresponding  crack 
surface  area  to  determine  the  quantity 

3 (6/P) 

3A 

needed  in  equation  (4.1.4)  for  calculation  of  the  energy 
release  rate.  A  sample  graph  for  compliance  against  crack 
surface  area  corresponding  to  a  specimen  of  thickness  3/8 
inch  is  shown  in  Figure  4.7.  Energy  release  rates  were  then 
calculated  using  equation  (4.1.4)  and  are  included  in  Table 
4.1. 

4.5  COMPARISON  OF  EXPERIMENTALLY  DETERMINED  ENERGY  RELEASE 
RATES  WITH  THE  HARTRANFT-SIH  MODEL 

Hartranft  and  Sih  [1]  carried  out  an  approximate  three- 
dimensional  analysis  for  a  finite  thickness  plate  containing 
a  through  crack  under  tensile  loading.  Figure  4.1.  In  this 
refined  theory,  the  rate  at  which  energy  is  released  with 

crack  extension  is  found  to  vary  through  the  plate  thickness. 
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An  average,  over  the  plate  thickness,  of  the  energy  release 
rate  per  unit  of  new  crack  lnegth  created  Is  given  In  the 
form 


G  =  f(h/2a,v,p,e)  7rg2a^1-~-v-2-^  (4.5.1) 

where  a  Is  the  applied  stress,  2a  is  the  crack  length,  h  Is 
the  plate  thickness,  v  Is  Poisson's  ratio,  and  E  Is  Young's 
modulus.  The  parameter  e  is  a  measure  of  the  boundary  layer- 
thickness  at  the  plate  surfaces  where  observed  behavior  Is 
qualitatively  from  that  of  the  specimen  interior.  These 
surface  layers  were  introduced  to  account  for  phenomena  such 
as  shear  lip  formation.  In  the  specimen  interior  the  varia¬ 
tion  of  the  near  tip  stress  field  with  the  thickness  coordi¬ 
nate  z  is  given  by  the  multiplicative  functional  coefficient 

f(z)  =  cospz 

The  value  of  p  is  to  be  chosen  from  experimental  observations. 

It  should  be  noted  that  if  the  crack  extends  without 
change  in  shape  of  its  leading  edge,  the  average  energy  re¬ 
lease  rate  per  unit  increase  of  crack  length  and  the  energy 
release  rate  per  unit  of  new  crack  surface  created  are  iden¬ 
tical.  With  this  assumption,  the  energy  release  rate  5 
(equation  (4.5*1))  and  the  experimentally  determined  energy 
release  rate  at  fracture  are  both  plotted  against  the  non- 
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dimensional  thickness  parameter  h/2a  in  Figure  4.9.  A  least 

squares  procedure  has  been  used  to  fit  the  experimental  curve 

through  the  data.  The  asymptotic  value  of  the  Hartranft-Sih 

h 

equation  for  G  at  large  ^  where  the  specimen  is  in  a  state 
of  plane  strain  and  S  ->  —  — I  has  been  chosen  to  agree 

with  the  experimental  results. 

In  two-dimensional  elasticity,  the  energy  release  rate  G 
is  used  as  a  failure  criterion  for  thick  specimens  where  plane 
strain  conditions  are  assumed  to  hold.  A  critical  value  of 
the  energy  release  rate  G^c  known  as  the  fracture  toughness 
is  obtained  experimentally  for  a  given  material.  Failure  of 
a  structural  component  of  the  same  material  containing  a  flaw 
is  predicted  at  the  level  of  loading  such  that 


G(load,  geometry)  ■+  G-^c 


Here,  an  attempt  is  made  to  extend  the  fracture  toughness 
concept  beyond  the  plane  strain  range.  A  thickness  dependent 
energy  release  rate  at  fracture  has  been  determined  experi¬ 
mentally  for  plexiglas.  Figure  4.9,  and  an  approximate  for¬ 
mula,  equation  (4.5.1),  can  be  used  to  determine  the  depen¬ 
dence  of  the  energy  release  rate  on  specimen  loading  and 
geometry.  The  G  obtained  in  equation  (4.5.1)  can  be  equated 
to  the  GeXp  for  given  h/a  to  determine  the  failure  load,  i.e.. 


P 

cr 


3exr9h2E 

-  r  exP — 

LafiT  ( 1-v^ )  J 


1/2 


(4.5.2) 


where  f  is  given  in  equation  (4.5-1).  These  calculations 
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have  been  carried  out  for  each  of  the  specimens  tested  and  a 
comparison  of  P  (exp)  with  P  (predicted)  is  given  in 

OX  ox 

Table  4.2. 

A  graphical  comparison  of  experimental  and  predicted 
values  for  P  is  shown  in  Figure  4.10.  The  quantity 

OX 

P  /a/h  is  plotted  against  the  thickness  parameter  h/2a. 

O  X 

This  parameter  (P  /a/h)  is  seen  to  be  relatively  insensitive 
to  h/2a  for  h/2a  >  .35  and  in  that  range  experimental  and  pre¬ 
dicted  values  are  in  good  agreement.  The  experimental  results 
for  the  range  h/2a  <  .35  are  observed  to  depend  strongly  on 
the  actual  plate  thickness  h  and  continuation  of  the  curve 
through  these  points  would  be  deceiving. 

In  considering  the  results  presented  in  Table  4.2  and 
Figure  4.10,  note  that  for  each  plate  thickness  tested  (with 
the  exception  of  h  =  0.5")  the  G  is  found  to  increase  with 
h/2a  until  a  maximum  value,  then  decrease  slightly  and  level 
out.  These  results  are  typical  of  those  obtained  many  ex¬ 
perimental  programs  in  fracture  mechanics.  The  level  portion 
of  the  curve  for  large  h/2a  is  known  as  the  plane  strain  re¬ 
gion.  The  peak  is  associated  with  plane  stress  and  the  left¬ 
most  portion  where  G  increases  with  h/2a  is  thought  to  be 
the  consequence  of  a  nonlinear  necking  phenomenon.  Tests 
corresponding  to  this  range  where  plasticity  is  considered 
to  be  important  and  the  Hartranft-Sih  elastic  model  is  not 
expected  to  be  accurate  are  marked  with  superscript  (*). 
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With  the  exception  of  the  results  for  0.5"  thickness 
specimens  which  are  suspect,  the  experimental  values  and  the 
predictions  for  critical  load  are  in  reasonable  agreement 
over  the  range  where  the  nonlinear  effects  are  thought  to 
be  negligible.  These  results  are  based  on  a  rather  arbitrary 
choice  for  the  parameters  p  and  e  contained  in  equation 
(4.5.1).  Further  theoretical  and  experimental  work  is  being 
carried  out  presently  to  develop  a  rational  procedure  for 
choosing  these  or  other  related  parameters.  Similar  sequences 
of  tests  are  needed  for  other  materials  to  substantiate  the 
present  results.  Even  with  these  limitations,  the  initial 
results  obtained  here  are  very  encouraging  in  that  they  sug¬ 
gest  a  possible  procedure  for  overcoming  the  testing  restric¬ 
tions  on  specimen  thickness  based  on  two-dimensional  plane 
strain  analyses,  and  for  using  fracture  mechanics  to  predict 
failure  of  three-dimensional  structural  components.  This 
kind  of  information  is  necessary  before  reasonable  strength 
analyses  can  be  carried  out  for  laminated  plates  containing 
flaws  where  the  stress  fields  near  the  flaws  are  expected  to 
be  three-dimensional  in  character. 
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TABLE  4.2  -  COMPARISON  OP  PREDICTED  AND  EXPERIMENTAL  VALUES 
FOR  FAILURE  LOAD 

G  P 

exp  cr 


h 

a 

lb-in/in2 

exp. 

pred. 

h/2a 

3/4  " 

3/4" 

7.45 

20 

( 10 ) 2lbs 

21.5(10) 2lbs 

1/2 

3/4  " 

5/8" 

7.0 

22.5 

(10) 2lbs 

22.5(10) 2lbs 

3/5 

3/4  " 

1/2" 

6.45 

28 

(10) 2lbs 

25  ( 10 ) 2lbs 

3/4 

3/4  " 

1/4" 

6.30 

34 

(10) 2 lbs 

35  ( 10 ) 2lbs 

3/2 

5/8  " 

3/4" 

6.50 

13.5 

(10) 2lbs 

16.4(10) 2lbs 

5/12 

5/8  " 

VJ1 

\ 

oo 

7.25 

18 

( 10 ) 2lbs 

19.6(10) 2lbs 

1/2 

5/8  " 

1/2" 

7-03 

24 

(10) 2lbs 

21.5(10) 2lbs 

5/8 

5/8  " 

1/4" 

6.80 

32 

(10) 2lbs 

29.2(10) 2 lbs 

5/4 

1/2  " 

3/4" 

7.80 

11 

( 10 ) 2lbs 

15.0(10) 2lbs 

1/3 

1/2  " 

5/8" 

7.40 

13.7 

(10) 2lbs 

16.0(10) 2lbs 

2/5 

1/2  " 

1/2" 

7.75 

22 

(10) 2lbs 

17.8(10) 2lbs 

1/2 

1/2  " 

1/4" 

6.30 

30 

(10) 2lbs 

23*0(10) 2 lbs 

1 

3/8  " 

3/4" 

8.60 

7-2 

(10) 2lbs 

11.5(10) 2lbs 

1/4 

3/8  " 

5/8" 

9-25 

9.5 

( 10 ) 2lbs 

12.0(10) 2lbs 

3/10 

3/8  " 

1/2" 

7.55 

11. 25 ( 10 ) 2lbs 

13.7(10) 2lbs 

3/8 

3/8  " 

3/8" 

7.0 

14.5 

(10) 2 lbs 

15.0(10)2lbs 

1/2 

3/8  " 

1/4" 

7.35 

18.7 

(10 ) 2lbs 

18.5(10) 2lbs 

3/4 

1/4  " 

3/4 

7.05 

4.55(10)2lbs 

8.2(10) 2  lbs 

1/6 

1/4  " 

1/2 

8.25 

10.8 

(10) 2lbs 

9. 6 (10) 2  lbs 

1/4 

1/4  " 

1/4" 

6.45 

13.8 

( 10 ) 2lbs 

13.2(10)2lbs 

1/2 

1/8  " 

3/4" 

6.97 

3.6 

( 10 ) 2lbs 

4.0(10) 2lbs 

1/12 

1/8  " 

5/8" 

7.25 

4.4 

(10) 2lbs 

4. 5 (10) 2  lbs 

1/10 

1/8  " 

1/2" 

6.40 

5-0 

(10) 2  lbs 

5. 0(10) 2  lbs 

1/8 

1/8  " 

3/8" 

6.65 

5-8 

(10) 2  lbs 

5.6(10 ) 2lbs 

1/6 

1/8  " 

1/4" 

6.72 

7.6 

(10) 2lbs 

7.0(10) 2lbs 

1/4 

1/16" 

3/4" 

5.80 

2.5 

( 10 ) 2lbs 

2 . 4 ( 10) 2lbs 

1/24 

1/16" 

5/8" 

6.20 

2.8 

(10 ) 2 lbs 

2. 6 (10) 2  lbs 

1/20 

1/16" 

1/2" 

7.02 

3.4 

(10) 2lbs 

2.9(10)2lbs 

1/16 

1/16" 

3/8" 

7.10 

4.5 

(10) 2lbs 

3.1(10)2lbs 

1/12 

1/16" 

1/4" 

5.75 

5.0 

( 10 ) 2lbs 

3.6(10)2lbs 

1/8 

SECTION  V 


A  PENNY-SHAPED  CRACK  IN  AN  ISOTROPIC  LAYER  SANDWICHED 
BETWEEN  TWO  ORTHOTROPIC  HALE-SPACES 

There  is  considerable  interest  in  the  possible  struc¬ 
tural  applications  of  fiber-reinforced  laminar  composites. 
These  composites  are  produced  by  adhering  sheets  of  fiber 
reinforced  material  to  create  a  laminar  structure.  Separa¬ 
tion  of  the  individual  layers  through  debonding,  a  process 
known  as  delamination,  is  perhaps  the  most  often  observed 
mode  of  failure  for  this  type  of  composite  system.  The  aim 
of  this  study  is  to  determine  the  influence  of  the  geometric 
and  material  properties  of  the  composite  on  the  loading  con¬ 
ditions  under  which  delamination  will  occur.  A  flaw  is  pos¬ 
tulated  to  exist  within  an  adhesive  layer  in  the  form  of  a 
penny-shaped  crack  parallel  to  the  plane  of  the  layers.  The 
shape  chosen  for  the  flaw  clearly  simplifies  the  mathematical 
analysis,  yet  it  also  models  the  physical  situation  realisti¬ 
cally.  Kassir  and  Sih  [32]  have  demonstrated  that  cracks  of 
arbitrary  shape  subjected  to  stress  fields,  which  are  approx¬ 
imately  uniform  relative  to  the  scale  of  the  flaw,  will  tend 
to  grow  into  the  shape  of  a  circular  disk  of  discontinuity. 

A  stress  analysis,  for  a  model  of  fiber-reinforced  composite, 
will  be  described  which  yields  a  measure  of  the  stress  inten¬ 
sification  along  the  edge  of  the  flaw  as  a  function  of  the 
thickness  of  the  bond  and  the  elastic  properties  of  both  the 
adhesive  and  the  fiber-reinforced  layers.  The  rate  at  which 

energy  could  be  made  available  for  the  creation  of  new  sur- 
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face  associated  with  enlargement  of  the  flaw  is  determined 
thereby  enabling  prediction  of  the  conditions  under  which 
delamination  will  occur  for  a  given  adhesive. 

5*1  MODEL  FOR  FIBER  REINFORCED  LAMINATE 

Consider  the  composite  made  up  of  a  number  of  fiber- 
reinforced  layers  bonded  together  by  relatively  thin  layers 
of  adhesive  materials.  One  of  the  adhesive  layers  contains  a 
penny-shaped  crack  which  is  parallel  to  the  bonding  surfaces. 
Figure  (5.1).  The  stated  purpose  of  the  analysis  is  to  de¬ 
scribe  accurately  the  solution  in  the  vicinity  of  the  flaw 
and  the  influence  of  the  geometric  and  material  properties  of 
the  composite  system  on^t-hat  solution.  With  this  in  mind, 
the  details  of  the  composite  structure  are  ignored  everywhere 
except  in  the  neighborhood  of  the  crack.  The  sequences  of 
reinforced  and  adhesive  layers  above  and  below  the  layer  which 
contains  the  flaw  are  idealized  as  two  identical,  transverse¬ 
ly  isotropic  half-spaces  with  preferred  directions  normal  to 
the  layers.  Figure  (5.2).  The  elastic  properties  of  these 
half-spaces  are  chosen  to  model  the  averaged  effects  of  the 
outer  layers  on  the  solution  in  the  vicinity  of  the  flaw. 

The  parameters  in  this  model  which  can  be  varied  to  represent 
the  influences  of  the  composite  structure  on  the  local  behav¬ 
ior  about  the  crack  are:  the  ratio  of  adhesive  layer  thick¬ 
ness  to  crack  radius,  the  Young's  modulus  and  Poisson's  ratio 
for  the  adhesive,  and  the  five  material  properties  of  the 

transversely  isotropic  half-spaces.  For  example,  the  relative 
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stiffness  of  the  half-spaces  parallel  to  the  bond  plane  as 
compared  to  their  stiffness  normal  to  that  plane  is  taken  to 
be  a  measure  of  the  stiffness  and  quantity  of  the  fiber  rein¬ 
forcements. 

Solutions  will  be  presented  for  the  problem  of  pressure 
loading  on  the  surfaces  of  the  penny-shaped  crack.  The  super¬ 
position  principle  for  linear  elasticity  can  be  used  to  demon¬ 
strate  that  the  sum  of  these  solutions  plus  the  corresponding 
ones  for  the  remote  tensile  loading  of  the  same  body  without 
the  crack  gives  the  solutions  for  remote  tensile  loading  nor¬ 
mal  to  the  bond  line  of  the  composite  containing  the  flaw. 

The  solution  for  the  remote  loaded  composite  which  does  not 
contain  a  crack  is  nonsingular  and  thus  will  not  effect  the 
total  solution  in  the  immediate  vicinity  of  the  crack  edge. 
Therefore,  as  far  as  the  near-tip  solution  is  concerned,  for¬ 
mulation  of  the  problem  for  either  pressure  loading  on  the 
crack  faces  or  remote  tensile  loading  lead  to  identical  re¬ 
sults.  The  former  is  amenable  to  solution  by  transform  pro¬ 
cedures  [33]  and  will  be  considered. 

5.2  FORMULATION  OF  GOVERNING  EQUATIONS 

For  convenience,  the  superscript  (1)  will  be  used  to  re¬ 
fer  to  the  isotropic  layer  containing  the  crack,  while  symbols 
with  a  superscript  (2)  will  be  associated  with  the  half-spaces. 
The  governing  equations  for  the  layer  and  half-spaces  cam  be 
expressed  in  terms  of  Airy  stress  functions  as 
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,92  ,19  ,  92  x  ,(1) 

(9rrr  +  r  3r  +  9z7‘)  *  =  0 


and  [34] 


(5.2.1) 


r9  +  I  1_  ,  1  3 2  \  / 3 2 

^ 3r2  r  3r  3z2  ' 3r' 


+  IL+1 

r  3r  S; 


92  U<2) 


where 


S,  =  a  +  3,  S0  =  a  -  $ 


/  M  a  _  y  P-l 

a  -  /  — ,  3  -  /  — 


p  _  1  rE^2  Vg^2^-2v^2^  (1+v^2^ )  , 

2  E(2)(1-v2(2)/E(2)) 

H  =  r  l-V(2)E<2>/j(2)  .1/2 

E<2)(l-v‘(2,)/E(2) 

(2 )  (2) 

Here,  E  and  v  are  the  Young's  modulus  and  the  Poisson's 
ratio,  respectively,  in  planes  parallel  to  the  interface, 
while,  E  ,  v  ,  and  Gv  are  the  Young's  modulus,  Poisson's 

ratio,  and  shear  modulus  in  the  transverse  direction  for  ma¬ 
terial  (2).  Solutions  will  be  presented  for  choices  of  ma¬ 
terial  properties  such  that  $  is  real. 


The  continuity  conditions  across  the  material  interface 


-97- 


^1}(h,r)  =  a<2)(h,r) 

T^}(h,r)  =  x£2)(h,r) 

(5.2.2) 

(h,r )  =  u£2)(h,r) 

(h,r)  =  uj; 2  ^  ( h , r ) 

The  mixed  boundary  conditions  on  the  z  =  0  plane  are 

(o,r)  =  -p,  r<a 

(o,r )  =  0,  r>a  (5-2.3) 

TP2 } (°>r )  =  0 

Let  G^^(z,£),  i  =  1,2,  be  the  Hankel  transforms  of  the 
Airy  stress  functions  <j>^^(z,r)  given  by 

G^(z,£)  =  /  r4>  ^ 1  ^  (z,r)JQ  (£,r)dr 
o 

Application  of  the  Hankel  transform  to  the  governing 
differential  equations  (5-2.1)  results  in  the  following  equa¬ 
tions  for  the  transforms  of  the  Airy  stress  functions: 

(gfr  -  £2)2G(1)(z,£)  =  0  (5-2.4) 
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%  &  -  «!)%  df*  -  52)Q(2)(Z,C)  =  0  (5.2.5) 

The  well-known  solutions  to  these  equations  are: 
G(1)(z,£)  =  (A(1)(£)  +  B(1) (£Kz)eosh£z 

+  (C(1)(?)  +  D(1) (C)?z)sinhCz  (5.2.6) 

G(2)(z,C  =  A(2)(£)e~?(a+0)z  +  B(2)(Oe"^(a_3)z 

+  C(2) (£)e^(a+3)z  +  D(2)  (Oe^a_B)z  (5.2.7) 


Here,  the  function  C^2^(£)  and  D^2^(£)  are  required  to  be 
identically  zero  so  that  the  solution  will  remain  bounded  for 
large  z. 

The  stress  components  and  the  displacements  for  the 
layer  and  the  half-space  can  be  expressed  in  terms  of  the 
transformed  stress  functions  G^^(z,?)  and  G^2^(^,z),  respec¬ 
tively  as 
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00 


^1}(z,r)  =  /  5CG(JL  -  ^TTT  S2G[d)]J0U,r)dS 


o  >zzz  1-v 


r^(z,r)  =  /  £2[  ~  7yy  G  zz 

rz  o  (l-vUJ)  jZZ 


2[  -  7TV  G(ll  +  ?2G(1)]J,(C,r)dC 


ui1}(z,r)  =  ,U+V  7TV  7  ?2G(d)J1(?,r)d? 


E^'CL-V^')  c 


uz1)(z>r)  =  (i)+V  TIT7  7  £[(i-2v(1))G(^ 

z  Eu;(l-Vu;)  0  ,ZZ 

-  2(l-v(1))C2G(1)]J0(C,r)d^ 


a^2)(z,r)  = 


-  /  E^zzz  - 


(5-2.8) 


+  ^  (b-l)?2G^}J1(C3r)dC 


a^2)(z,r)  =  7  UdG(2L  -  C?2G(2)  ]Jn(?,r)d£ 


j z  J"0 


Trz’<z>r)  = 


-  /  5a[aO(^  -  52Q(2)]J1(5Jr)dC 

O  5 


u<2)(z,r)  =  — 1  / 

E  o 


/...  rN  7  ^r(2v  2  a+d)  d(2)  r  ( 1+b  )  V  2  +C  -|  „ 
(z’r>  -  J  -  G  -  L  -T^n J  x 


X  S2G(2)}J0U,r)d£ 
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where 


v^2^(l+v^2^) 

=  '  i-ot2)E<2Vi(2) 

b  =  [v<2)  -  v(2)(e(2)/b(2)  -  v<2>E(2)/i(2))]/  x 
X  [1  -  v2 (2 )e^2 )/E^2  ^ ] 

c  =  [E(2)/G(2)  -  v(2)(l+v(2))]/[l  -  V2(2)E(2)/E(2)] 

d  =  [E(2)(l  -  v2(2))/E(2)]/[1  -  v2(2)E(2)/E(2)] 

The  general  solutions  for  G^(z,5)>  equation  (5.2.6), 
and  equation  (5*2.7)  are  substituted  into  the  above  expres¬ 
sions  for  the  stress  and  displacement  components  and  the 
boundary  conditions,  equations  (5.2.2)  and  (5.2.3)  are  ap¬ 
plied.  Five  of  the  unknown  functions  of  £  can  be  expressed 
in  terms  of  the  sixth,  reducing  the  problem  to  the  set  of 
dual  integral  equations  with  the  mixed  boundary  conditions  on 
the  z  =  0  plane. 

/  £3A(1)U)J 0(5,r)d£  =  0,  r>a 

o 

(5.2.9) 

/  S4f  (£)A^(£  )JQ  (£,r)d£  =  -p,  r<a 

o 

where  f(£),  a  complicated  function  of  the  geometric  and  ma¬ 
terial  parameters  of  the  composite,  is  given  in  the  Appendix 
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B.  The  other  unknown  functions  B^^(£)>  D^^(£)j 

A^^),  B^^(£)  can  be  expressed  in  terms  of  A^\e,)  ,  see 
the  Appendix  B. 

Solution  procedures  have  been  developed  for  this  set  of 
dual  integral  equations  [28]  and  the  method  of  solution  will 
be  briefly  outlined  here  for  completeness.  The  first  of  equa¬ 
tions  (5.2.9)  is  satisfied  by  introducing  the  function 

W(r)  =  /  £3A(1)(S)J0(5,r)dS  (5.2.10) 

o 

with  the  condition  W(r)  =  0  for  r>a.  The  Hankel  transform 
is  applied  to  equation  (5.2.10)  to  obtain 

52A^(0  =  /  rW(r)JQ  (£  ,r)dr 
o 

The  second  equation  (5-2.9)  may  thus  be  written  as 

00 

/  £2f (0{nW(n)J0(£n)dn}Jo(£>r)d£  =  p(r)>  r<a  (5. 2. 11) 

o 

Since  W(x)  is  associated  with  the  displacement  field  it  must 

1/0 

be  bounded  at  x±a  and  of  order  (x±a)  '  [28]  as  |x|+a,  i.e., 

W(x)  =  (j)  (a)  (a2-x2 ) 1//2  +  ... 

For  this  reason,  W(r)  is  expressed  in  the  form 
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/  ~~~~dT ,  0<r<£ 


W(r )  = 


,  r>a 


where  4>(t)  is  continuous  on  the  interval  [0,a]. 


The  expression  for  W(r)  is  substituted  into  equation 
(5*2.11).  The  resulting  equation  is  multiplied  by  r  and  then 
integrated  in  r  from  r  =  0  to  r  =  R.  After  a  formal  change 
in  order  of  integration  and  application  of  the  identity 


x  CnJn(CJn)dn 

/  - — — ■■■■ — dri  =  sin(£t) 


the  following  equation  is  obtained 


a  00  R 

R  /  t <J> ( x )  /  f(?)sin(?T)J1(?R)dCdT  =  /  p (r )rdr 
0  0  o 


where  J(£R)  is  the  first  order  Bessel  function  of  the  first 
kind. 


For  the  special  case  of  constant  pressure  the  above 
equation  is  reduced  to: 


/  /  f  (C)sin(£x)J,  (£R)d£dT  = 


(5.2.12) 


By  making  the  substitution 


f*(£)  -  -  f(S)2v(1) (l-v(1) )  -  1.0 
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and  applying  the  identity 


-104- 


t  =  at,  n  =  as,  t<J>(t)  =  -  p’a$(x) 


lead  to 


$(t)  +  /  $(s)k(s,t)ds  =  t 
0 


(5.2.14) 


with  the  symmetric  kernel 


:(s,t)  =  f  /  f*(|)sin(£s)sin(£t)d£,  0<t<lj  0<s<l 


The  equation  (5-2.14)  is  solved  numerically  for  $(t). 

The  function  $(t)  can  be  determined  exactly  for  the  case 
of  a  limitingly  thin  adhesive  layer.  Note  that  the  function 
f(5)  reduces  to  a  constant,  say  f  ,  for  h  =  0.  Refer  back 
to  equation  (5-2.12),  divide  by  f  and  use  the  identity, 
equation  (5.2.13),  to  reduce  equation  (5.2.12)  to  the  Abel 
integral  equation 


/  ,  E«i  r<6 


Solution  of  this  equation  and  application  of  the  normaliza¬ 
tion  described  above  gives 


$(t)  =  -  — t—c - 7^-r - ,  for  h  =  0 

2vu;(l-vu;)f 
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Returning  to  the  solution  of  the  original  set  of  dual 
integral  equations,  A^^(E)  takes  the  form 

A(1)(?)  =  -  ®(l)cos(5a) 

1  i  a 

+  /  tt  cos(^at)^r-  [$(t)]dt  (5.2.15) 

0  ^ 

where  $(1)  is  to  be  interpreted  as  the  limit  of  $(t)  as  t 
approaches  1.0  from  below.  Once  the  function  $(t)  is  evalu¬ 
ated  from  the  Fredholm  integral  equation,  equation  (5-2.15) 
may  be  substituted  into  equations  (5-2.6)  and  (5-2.7)  to  ob¬ 
tain  and  from  which  the  displacements  and  stresses 

at  every  point  may  be  calculated. 

5-3  LOCAL  STRESS  FIELD 

It  is  not  difficult  to  show  that  the  integral  expression 
in  equation  (5.2.15)  remains  bounded  at  the  crack  edge  r  =  a 
and  that  the  contribution  to  the  singular  behavior  of  the 
stresses  is  fully  described  by  the  leading  term  in  equation 
(5.2.15).  The  desired  explicit  representation  of  the  portion 
of  the  solution  that  becomes  unbounded  at  the  penny-shaped 
crack  border  may  be  obtained  by  recourse  to  Laurent  series 
expansions  and  the  use  of  elementary  Bessel-integral  identi¬ 
ties.  Adopting  polar  coordinates  P  and  6,  as  shown  in  Figure 
(5- 3) ,  it  is  found  that 
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ar  =  k1(2p)-1,/2  cos  [1  -  sin  j  sin  |-0] 
az  =  k1(2p)-1//2  cos  |-  [1  +  sin  J-  sin  |-0] 

Trz  =  k1(2P)_1//2  cos  §  ^sin  f  cos 
where  k^  is  the  stress  intensity  factor  defined  as 

kl  =  f  $ (l)p/a 
5.4  NUMERICAL  RESULTS 

Numerical  solutions  for  the  stress  intensity  factor 
were  obtained  for  the  special  case  =  v^2^  =  v^2^  =  3, 

G^2  ^  =  E^ 2  V2  (1+v  ^2  ^ ) .  For  this  choice  of  G^2^,  3,  equation 

rp(2)  ~  /  p  \  2 

(5-2.1)  is  real  for  1  <  ^ ^  <  [l/vv  ']  and  the  general 

solution  to  equation  (5-2.5)  for  the  transformed  Airy  stress 

function,  equation  (5-2.7)  is  correct.  Alternate  expressions 
~  ( 2 ) 

for  G  are  possible  which  results  in  real  values  for  3  over 
other  ranges  of  E'  '/Ev  . 

Numerical  results  for  the  variation  of  stress  intensity 
factor  over  the  range  of  material's  properties  discussed 
above  are  exhibited  in  Figures  (5.4)  to  (5.6).  In  these 
plots,  the  stress  intensity  factor  has  been  normalized  by  its 
value  for  a  penny-shaped  crack  of  the  same  radius  and  sub¬ 
jected  to  the  same  pressure  loading  but  contained  in  an  in- 
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finite  isotropic  body  (k1  =  p/a)  [24].  Thus  values  of  the 
normalized  stress  intensity  factor  which  are  less  than  1.0 
indicate  that  less  load  is  transferred  to  the  crack  tip  by 
the  laminar  composite  than  would  be  transferred  in  an  equiva¬ 
lent  homogeneous,  isotropic  body.  The  dependence  of  the 
stress  intensity  factor  on  the  degree  of  anisotropy  of  the 
center  layers  E^/E^  is  illustrated  in  Figures  (5-4)  and  (5-5) 
for  various  values  of  the  relative  stiffness  of  the  cracked 
layer  to  that  of  the  half  spaces  (in  the  direction  normal  to 
the  crack  face),  E-JE^  at  a  number  of  distinct  ratios  of 
layer  thickness  to  crack  diameter.  Increasing  the  relative 
stiffness  of  the  half  spaces  parallel  to  the  crack  face  re¬ 
duces  the  load  transfer  to  the  crack  tip.  Further,  any 
deviation  of  the  stress  intensity  factor  from  the  one-mate- 
rial  result  (k^  =  —  p/a)  is  amplified  as  the  layer  thick¬ 

ness  to  crack  diameter  ratio  h/a  is  decreased.  This  effect 
is  demonstrated  in  Figure  (5*6).  It  should  be  noted  that  the 
slope  of  the  stress  intensity  factor  layer  thickness  curves 
is  infinite  at  h/a  =  0.  This  mathematical  result  has  been 
demonstrated  earlier  in  connection  with  a  different  physical 
problem  [28]. 

5.5  DISCUSSION 

The  primary  purpose  of  these  calculations  is  associated 

with  the  use  of  this  boundary  value  problem  as  a  model  to 

study  the  effect  of  the  geometric  and  material  properties  of 

the  fiber  reinforced,  laminar  composite  system  on  a  flaw  in 
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an  adhesive  bond.  With  this  interpretation  in  mind,  the  pre¬ 
dictions  are  that  load  transfer  to  the  flaw  is  decreased  as 
the  relative  stiffness  of  the  layer  material  to  that  of  the 
adhesive  is  increased  and,  in  particular,  as  the  volume  per¬ 
centage  or  stiffness  of  the  reinforcing  fiber  is  increased, 
producing  a  higher  apparent  stiffness  for  the  layers  parallel 
to  the  interfaces.  The  effects  are  more  pronounced  when  the 
adhesive  layers  are  thin  or  the  flaw  is  large. 


These  estimates  of  load  transfer  can  be  related  to  the 
prediction  of  conditions  for  delamination  via  a  knowledge  of 
the  energy  necessary  to  create  new  surface  for  the  bonding 
material.  The  rate  at  which  energy  would  be  made  available 
for  the  creation  of  new  surface  if  the  crack  expanded  uniform¬ 
ly  is  readily  attainable  from  the  numerical  results  presented 
above  via  the  energy  release  rate  equation  for  cracks  in  an 
isotropic  material  [24]: 


G 


4(l-v2(1)) 


k 


2 

1 


Thus  crack  extension  is  predicted  when  G  is  greater  than  the 
corresponding  material  property  Gc  =  2iTa2Y  of  the  adhesive. 
This  analysis  thus  provides  a  quantitative  criterion  for 
choice  of  bonding  material  in  terms  of  its  stiffness  and  en¬ 
ergy  required  to  create  new  surface. 
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APPENDIX  A  -  PARAMETERS  ARISING  IN  SOLUTION  FOR  IN-PLANE 
STRESS  VARIATIONS 


a4 

kl  =  2al 


an  ar 

■  lk3  =  ^  (2a^ 


'  2  *6  [(k^'+k* 

2  al  k3<i6+al> 

2 


Cg  -  1^^  -  2k^a g 


I,  i 

c3  =  ^  (k3b2-k1b1)  -  i  (k*-k*)a6 


d1  =  2k3ag  - 


d2  =  k1k3ag  +  ^  (k^^+k^bg) 


I3[^  k1k3b1  +  |  (k*+k*)b2] 


Tj-  k3(k2-3k2)a6 
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APPENDIX  B  -  FUNCTION  f(g)  RELATED  TO  THE  SECOND  EQUATION 
(5.2.9) 


f(C)  = 


1 _  H(g) 


where : 
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Figure  3- 3- 3-1  -  Dimensionless  Stress  Intensity  Factor  as  a  Function  of  Normalized 

Layer  Thickness  for  Formulation  without  Boundary  Layer 


Figure  3- 3 -3. 3  -  Dimensionless  Stress  Intensity  Factor  as  a  Function  of  Normalized 

Layer  Thickness  for  Formulation  without  Boundary  Layer 
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Figure  3 . 3 • 3 • 5  -  Dimensionless  Stress  Intensity  Factor  at  Various  Ratios  of  Material 

Stiffness  for  Formulation  including  Boundary  Layers 


Figure  3  -  3  -  3  -  7  -  Dimensionless  Stress  Intensity  Factor  at  Various  Ratios  of  Material 

Stiffness  for  Formulation  including  Boundary  Layers 
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Figure  3. 3. 3- 8  -  Dimensionless  Stress  Intensity  Factor  at  Various  Ratios  of  Material 

Stiffness  for  Formulation  including  Boundary  Layers 
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Figure  3. 3* 3- 9  -  Dimensionless  Stress  Intensity  Factor  at  Various  Ratios  of  Material 

Stiffness  for  Formulation  including  Boundary  Layers 
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Figure  3.3-3.11  -  Influence  of  Relative  Layer  Thickness  on  Normalized 

Stress  Intensity  Factor 
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Figure  3 .3.3 .12  -  Influence  of  Relative  Layer  Thickness  on  Normalized 

Stress  Intensity  Factor 


Figure  3.3.3* 13  -  Influence  of  Relative  Layer  Thickness  on  Normalized 

Stress  Intensity  Factor 
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Figure  4.2a  -  Plate  Containing  a  Central  Crack 

Figure  4.2b  -  Load-Displacement  Curve  for  Constant  Displacement 
During  Crack  Extension 

Figure  4.2c  -  Load-Displacement  Curve  for  Constant  Load  During 
Crack  Extension 


Figure  4.4  -  Testing  Details 
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Figure  4.6  -  Load-Displacement  Curves 


II 


5/  «  6 /  » 

:  732  732 


12/  " 

732 


15/  » 
o2 


Figure  4.8  -  Compliance  Curves  for  3/8"  Specimens 
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Figure  5*2  -  Remote  Loading  for  the  Composite  Cylinder 
Containing  a  Penny-Shaped  Crack  in  the 
Adhesive  Layer 
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Analyses  for  laminar  plates  are  reviewed  in  this  report.  It  is  observed  that 
while  approximate  global  stress  analyses  have  been  performed  for  laminates,  little  has 
been  accomplished  in  the  area  of  strength  analysis.  The  work  presented  here  is  con¬ 
centrated  in  that  area,  i.e.,  mathematical  models  to  consider  the  influence  of  flaws 
on  failure  loads  associated  with  crack  propagation  and/or  delamination.  In  particular, 
an  approximate  three-dimensional  formulation  for  laminates  based  on  the  variational 
principle  of  minimum  complementary  potential  energy  is  presented  and  applied  to  the 
problem  of  a  three  layer  composite  plate  containing  a  through  crack.  This  model  is 
employed  to  study  the  influences  of  the  geometric  and  material  characteristics  of  the 
composite  on  the  local  stress  field  in  the  vicinity  of  the  flaw. 

Results  of  an  experimental  program  to  determine  the  influence  of  plate  thickness 
on  critical  load  for  a  centrally  cracked  tensile  specimen  are  reported.  Extension  of 
the  experimental  program  to  consideration  of  laminated  specimens  is  currently  in  prog¬ 
ress. 

A  stress  analysis  has  been  carried  out  for  another  problem  of  particular  interest 
in  the  study  of  failure  for  composite  systems.  The  effect  of  voids  in  adhesive  layers 
on  the  behavior  of  multi-layer  fibrous  composites  is  modeled  by  the  problem  of  a  penny 
shaped  flaw  contained  in  an  isotropic  layer  joining  two  orthotropic  half  spaces. 
Stresses  in  the  vicinity  of  the  crack  edge  are  calculated  in  terms  of  the  loading  and 
the  properties  of  the  components  for  this  composite  system. 
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